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I.  INTRODUCTION 


With  the  discovery  of  the  infrared  radiation  by  Herschel,  its  applications  became  a 
major  technology.  Since  then,  the  designers  and  the  engineers  have  searched  for  a 
summary  measure  of  quality  in  the  evaluation  process  of  infrared  systems.  Finally,  the 
concept  of  minimum  resolvable  temperature  difference  (MRTD)  proposed  by  Sendall  and 
Lloyd  [Ref.  1]  has  been  widely  accepted  as  a  figure  of  merit  in  predicting  thermal  imaging 
system  (TIS)  performance  because  it  is  the  first  comprehensive  measure  that  includes 
system  resolution,  noise,  and  the  observer  into  calculations. 

The  first  standard  model  was  the  U.S.  Army  Night  Vision  Laboratory  (NVL) 
Static  Perfonnance  Model  (also  known  as  the  Ratches  Model),  which  was  developed  in 
1975.  It  was  a  one-dimensional  model  that  was  used  for  the  first-generation  serial  and 
parallel  scanning  thermal  imaging  systems.  The  model  gained  widespread  acceptance  for 
successfully  treating  the  human  eye-brain  recognition  process  as  a  matched  filter; 
however,  it  did  not  incorporate  sampling  effects  and  noise  sources  other  than  the  detector 
noise  [Ref.  1].  Moreover,  improvements  in  the  second-generation  thennal  imaging 
systems,  such  as  the  widespread  use  of  focal  plane  arrays  (FPA)  and  digital  signal 
processing  technology  made  it  necessary  to  revise  the  so-called  Ratches  Model. 

Therefore,  the  U.S.  Anny  Night  Vision  Laboratory  and  Electronic  Sensors 
Directorate  (NVESD)  presented  a  new  model,  called  FLIR-90,  that  was  later  refined  to 
FLIR-92.  FLIR-92  was  capable  of  introducing  the  3-D  noise  effects  into  the  model  and 
predicting  the  two-dimensional  MRTD  values  successfully,  but  it  did  not  account  for 
scene -phasing  effects.  Furthermore,  it  did  not  address  the  aliasing  problem  that  was 
experienced  due  to  the  undersampling  process  in  the  FPAs.  Since  the  treatment  of 
sampling  effects  limited  the  MRTD  predictions  to  sub-Nyquist  frequencies,  the 
information  above  the  Nyquist  frequency  was  not  completely  quantified  by  FLIR-92. 
Hence,  the  range  performance  of  the  system  was  artificially  limited  [Ref.  2]. 
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The  above  named  models  named  above  made  some  assumptions  regarding  human 
eye-brain  judgement  in  the  recognition  process,  which  brought  some  limitations  to  the 
performance  predictions  of  the  thermal  imagers.  Therefore  a  visibility  model,  which 
predicts  the  MRTD  of  second  generation  thennal  imaging  systems  based  on  a  minimum 
threshold  input  contrast  parameter  and  a  contrast  reduction  factor  due  to  aliasing  and 
blurring  effects  was  proposed  by  A.W.Cooper  and  R.J.Pieper  of  the  Naval  Postgraduate 
School  in  1994.  The  model  makes  no  assumptions  about  the  recognition  process,  and  thus 
it  suggests  defining  an  individual  threshold  value  for  each  observer  (either  human  or 
machine).  Moreover,  the  model  includes  the  aliasing  phenomenon  and  provides  MRTD 
predictions  beyond  the  Nyquist  frequency  limit  [Ref.  3]. 

In  1999,  NVTherm  was  developed  by  NVEOD  as  a  new  model  adopting  two 
major  modifications  to  FLIR-92,  which  were  the  new  eye  model  and  the  MTF  squeeze 
phenomenon.  The  new  eye  model  brought  a  limitation  to  the  eye  spatial  resolution 
capability,  which  made  the  results  more  realistic.  On  the  other  hand,  the  MTF  squeeze 
phenomenon  degraded  the  system  MTF  to  account  for  aliasing  effects  by  taking  the 
frequencies  beyond  the  Nyquist  limit  into  consideration  [Refs.  4  and  5]. 

There  is  no  doubt  that  staring  systems  will  be  dominant  in  this  discipline 
someday.  Furthermore,  many  possible  developments  are  expected  in  the  next  generation 
of  thermal-imaging  technology.  Today  some  of  the  topics  of  interest  are:  a  near-unity  fill- 
factor  detector  array,  two  or  more  colored-images,  automated  motion  detection,  tracking 
and  identification  functions  on  the  focal  plane  arrays,  user-friendly  image  processing,  3- 
D  image  processing  and  passive  millimeter  wave  imaging  on  the  same  FPA  with  infrared. 
[Ref.  6]  No  matter  how  far  the  thermal  imaging  system  technology  improves,  the 
practical  modeling  efforts  in  predicting  its  perfonnance  will  definitely  continue. 

The  main  objective  of  this  thesis  is  to  develop  a  virtual  thermal  image-processing 
model  [Ref.  23],  which  takes  some  of  the  current  thermal  imaging  system  modeling 
concerns  into  account,  such  as  noise,  sampling  and  aliasing,  and  creates  visual  images 
that  can  be  obtained  with  the  actual  thennal  imaging  systems  being  modeled.  A  further 
analysis  of  these  images  will  provide  some  useful  insights  into  the  problems  encountered 
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in  modeling  efforts.  Besides  that,  it  is  also  intended  to  validate  the  performance  of  the 
Visibility  Model  by  comparing  the  subjective  and  objective  MRTD  predictions  to  the 
measurement  results  in  this  thesis. 

As  a  starting  point,  the  fundamentals  of  thermal  imaging  systems  and  their 
operation  will  be  explained  in  Chapter  II.  Chapter  III  will  provide  a  general  overview  of 
the  principles  behind  the  FLIR-92  model.  Considering  that  the  Virtual  Thennal  Image- 
Processing  model  is  based  on  the  Visibility  Model,  Chapter  IV  will  discuss  the  Visibility 
concept  with  the  latest  modifications.  Then  the  experimental  setup  and  the  methodology 
used  in  subjective  and  objective  MRTD  measurements  will  be  presented  in  Chapter  V. 
This  chapter  will  also  provide  the  measurement  results  and  a  comparative  analysis  to  the 
predictions  obtained  from  the  Visibility  and  FLIR-92  Models.  Next,  the  Virtual  Thennal 
Image-Processing  model  will  be  introduced  in  Chapter  VI.  Chapter  VII  will  analyze  the 
images  created  by  this  model,  and  this  work  will  conclude  with  some  conclusions  and 
suggested  directions  for  future  studies. 
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II.  THERMAL  IMAGING  SYSTEMS 


The  purpose  of  this  chapter  is  to  give  a  summary  of  the  basic  principles  of  thermal 
imaging  system  operation,  which  will  be  necessary  for  the  reader  in  understanding  the 
content  of  the  following  chapters. 

A.  TYPICAL  THERMAL  IMAGING  SYSTEM  OPERATION 

The  region  of  the  electromagnetic  spectrum  covering  optical  wavelengths  is  called 
the  optical  spectrum,  including  the  ultraviolet  (UV),  the  visible,  and  the  infrared  (IR) 
spectral  bands.  The  ultraviolet  (UV)  spectral  band  occupies  the  wavelengths  from  0.1  to 
0.4  pm  while  the  visible  spectral  band  occupies  the  wavelengths  from  0.4  to  0.7  pm 
where  the  human  eye  can  produce  images.  Bounded  by  the  visible  and  the  microwave 
spectral  bands,  the  infrared  spectral  band  is  also  divided  into  four  sub-regions.  These  are: 
the  near  infrared  or  short- wavelength  infrared  (SWIR)  band  from  0.7  to  3  pm,  the  mid¬ 
wavelength  infrared  (MWIR)  band  from  3  to  8  pm,  the  long-wavelength  infrared  (LWIR) 
band  from  8  to  14  pm,  and  the  far  and  the  extreme  infrared  bands  from  14  to  1000  pm, 
respectively.  The  corresponding  spectra  can  be  seen  in  Figure  2.1  [Ref.  7]. 
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Figure  2.1.  The  Electromagnetic  Spectrum  [Ref.  8]. 
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Thermal  imagers  are  passive  optical  systems,  which  are  used  to  extend  the  human 
vision  beyond  the  visible  spectral  region  limits  by  making  the  naturally  generated 
radiation  in  the  mid-wavelength  and  long-wavelength  infrared  regions  visible.  Thennal 
imagers  (i.e.  Forward  Looking  Infrared  (FLIR)  Systems)  sense  the  self-emitted  thermal 
energy  of  the  target  and  the  background  in  the  field  of  interest.  Unlike  thermal  imagers, 
the  electro-optical  systems  (e.g.  night  vision  devices)  utilize  the  reflected  energy  from  the 
target  and  the  background  in  the  near  infrared  region,  which  is  generated  by  an  external 
energy  source  (i.e.  starlight).  Therefore,  the  reflectivity  characteristics  of  the  target  and 
background  become  an  interest  area  of  electro-optical  imaging  studies,  while  the 
emissivity  reflectivity  characteristics  of  the  target  and  background  become  an  interest 
area  of  the  thermal  imaging  studies  [Refs.  9  and  10]. 


Target  and  background 


Figure  2.2.  Thermal  Imaging  System  Block  Diagram  [Ref.  9]. 


A  simplified  diagram  illustrating  the  fundamental  components  of  a  thermal 
imaging  system  is  presented  above  in  Figure  2.2.  First,  a  target  should  produce  a 
sufficient  target-to-background  apparent  temperature  difference  in  order  to  be  detected, 
recognized  or  identified.  Then  the  radiation  from  the  scene  of  interest,  which  is  a 
combination  of  the  target,  the  background  and  the  atmospheric  properties,  passes  through 
the  atmosphere  on  the  way  to  the  sensor.  The  atmosphere  modifies  the  amplitude  and  the 
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phase  of  the  radiated  light  between  the  scene  and  the  sensor  by  absorbing  and  scattering. 
These  phenomena  together  degrade  the  information  content  of  the  radiation  reaching  the 
sensor.  Moreover,  unwanted  rays  of  light  are  reflected  into  the  sensor  aperture.  [Refs.  9, 
10,11] 

The  first  component  of  the  sensor,  the  optics,  gathers  the  incoming  radiation, 
filters  it  spectrally  and  spatially  and  focuses  it  at  the  optics  image  plane,  where  a  detector 
is  mounted.  The  detector  converts  the  optical  signals  into  analog  electrical  signals,  which 
represent  the  spatial  distribution  of  the  radiant  flux  intensity  leaving  the  scene.  [Ref.  10] 

The  sensor  detector  may  be  a  single  detector,  a  linear  array  (a  column  or  a  few 
columns)  or  a  two-dimensional  staring  array.  No  matter  what  kind  of  detector  system  is 
mounted  at  the  image  plane,  the  sampling  phenomenon,  which  is  necessary  for  the  image 
formation  in  the  latter  stages  of  the  process,  is  accomplished  there.  In  the  scanning 
systems  (single  detector  systems  or  linear  array  detector  systems)  the  image  is  moved 
across  the  detector(s)  while  the  image  is  sampled  perpendicular  to  the  scanning  direction. 
The  output  of  every  single  detector  symbolizes  the  local  scene  intensity  across  a  scanning 
line.  The  staring  array  technology  does  not  require  a  scanner  mechanism  because  the 
image  of  the  scene  is  sampled  by  the  physically  separated  detector  elements  in  the  array 
[Refs.  4  and  9]. 

Due  to  the  differing  objects  in  the  scene  of  interest,  the  scene  intensity  distribution 
results  in  various  detector  output  voltage  levels.  Therefore,  the  output  voltage  distribution 
of  the  detector  array  produces  an  electronic  signal  after  being  processed  and  amplified  by 
the  electronics.  The  electronic  analog  signal  is  usually  digitized  for  ease  of  manipulation 
in  the  reconstruction  process  of  the  image.  The  main  goal  of  the  electronics  here  is  to 
arrange  the  signal  at  the  detector  output  so  that  it  becomes  applicable  to  be  displayed  for 
the  human  observer.  The  video  display  may  or  may  not  be  a  component  of  the  system  in 
case  the  signal  is  arranged  for  an  automated  image  processing  system  (i.e.  automated 
target  recognizer  (ATR))  whose  output  is  adopted  for  operator  tracking  [Refs.  4,  9  and 
11]. 
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B.  LINEAR  SYSTEM  THEORY 

A  system  is  the  mapping  of  an  input  function  into  an  output  function.  The  input 
and  output  functions  have  different  characteristics  for  the  various  systems.  For  example, 
these  functions  for  the  electrical  networks  are  real  for  which  time  is  the  only  independent 
variable.  On  the  other  hand,  for  the  imaging  systems  where  space  is  the  two-dimensional 
independent  variable,  they  can  be  either  real  or  complex  [Ref.  12] 

A  convenient  mathematical  representation  of  an  imaging  system  can  be  the 
operator  S{}.  The  response  of  the  system  to  the  input  function  yields  the  output  function. 
If  i(x,y)  and  o(x,y)  are  the  input  and  the  output  functions  respectively,  the  input-output 
relation  of  the  system  can  be  represented  as: 

S(i(x,y)}  =  o(x,y)  (2.i) 

where  x  and  y  are  the  spatial  variables  [Ref.  9]. 

In  order  to  call  a  system  linear  there  are  two  essential  requirements  to  be  met:  the 
superposition  and  scaling  properties  should  be  applicable  to  the  system.  The 
superposition  property  suggests  that  the  response  of  the  system  to  the  sum  of  all  input 
functions  is  equal  to  the  sum  of  the  individual  system  responses  to  each  input  function  as 
stated  in  Equation  (2.2)  where  ii(x,y)  and  i2(x,y)  are  individual  input  functions. 

s {ij(x,  y)  +  i2 (x,  y )}  =  S {ix (x,  y )}  +  S{  i2 (x,  y)}  (2.2) 

On  the  other  hand,  scaling  property  suggests  that  a  scale  factor  applied  at  the 
input  results  in  the  same  scale  factor  at  the  output  of  the  system.  Using  the  same  notation, 
this  property  can  be  expressed  as  in  Equation  (2.3)  where  K  represents  the  constant 
scaling  factor  [Refs.  9  and  12]. 

S {K  i(x,  y)}  =  K  S {i(x,  y)}  =  o(x,  y)  (2.3) 
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The  advantage  of  linearity  is  essential  in  the  imaging  system  analysis  because  it 
allows  us  to  decompose  a  complicated  input  function  into  a  group  of  much  simpler  input 
functions  to  calculate  the  system  responses  to  each  of  these  elementary  inputs.  Thus,  the 
superposition  of  the  calculated  individual  responses  yields  the  total  system  response.  One 
such  decomposition  process  can  be  performed  using  Fourier  analysis  [Ref.  12]. 

The  system  response  to  a  delta  function  is  called  the  impulse  response  for 
electronic  networks,  or  point  spread  function  (PSF)  for  the  optical  imaging  systems. 
Since  imaging  is  a  convolutionary  process,  once  the  imaging  system’s  response  behavior 
to  a  delta  function  is  fully  described,  then  the  output  to  any  input  function  can  easily  be 
determined.  Equation  2.4  represents  the  two-dimensional  convolution  integral  of  the 
input  function  with  the  system’s  point  spread  function  h(x,y)  where  a  and  p  are  dummy 
variables. 

+a+a 

o(x,  y)  =  J  J  i(<r,  J3)  h(x  -  a,  y  -  j8)  dadjS  (2.4) 


The  Fourier  Transform  of  the  point-spread  function  is  called  the  optical  transfer 
function  (OTF)  while  the  Fourier  Transform  of  the  impulse  response  is  called  the  transfer 
function  in  a  general  fashion.  Note  that  any  transfer  function  is  the  frequency  domain 
representation  of  an  impulse  response  function.  Remembering  the  Fourier  property  that 
the  convolution  in  the  temporal  domain  corresponds  to  multiplication  in  the  frequency 
domain,  Equation  2.4  can  be  rewritten  as  Equation  2.5  where  0(c,p)  and  I(c,q)  are  the 
output  and  the  input  spectrum  representations  respectively,  and  H(£,r|)  is  the  transfer 
function.  Notice  that  I(^,p)  and  0(c,r|)  are  the  Fourier  Transforms  of  i(x,y)  and  o(x,y) 
functions  where  \  and  q  are  horizontal  and  vertical  spatial  frequencies  respectively  which 
have  units  in  cycles  per  radian. 

0(%,ri)  =  (2.5) 

At  this  point,  it  is  appropriate  to  discuss  the  relationship  between  temporal  and 
spatial  frequencies.  Just  as  temporal  frequency  is  often  a  convenient  parameter  in 
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describing  electrical  signals,  spatial  frequency  is  often  used  to  describe  objects,  their 
images  and  imaging  systems  designed  to  render  these  images.  Whereas  an  electrical 
signal  can  be  described  in  tenns  of  its  amplitude,  characteristics  or  in  the  Fourier 
Transform  domain  by  its  frequency  content,  an  object  can  be  described  in  space  by  its 
radiance  or  its  spatial  frequency  content. 

As  mentioned  previously  in  Chapter  II,  the  detector  component  converts  the 
incoming  photon  energy  into  electrical  energy.  This  suggests  the  necessity  of  using  both 
spatial  and  temporal  frequencies  in  thermal  imaging  system  analysis.  However,  spatial 
frequency  analysis  is  easier  to  apply  in  sensor  perfonnance  calculations  because  spatial 
frequency  is  range-independent  by  definition.  Therefore,  the  temporal  frequencies  are 
generally  converted  into  spatial  frequencies  using  Equation  2.6  for  scanning  systems. 

f  =  ZXVSCAN  (2.6) 


Where 

/  is  the  temporal  frequency  in  Hz 

vscan  is  the  scan  velocity  in  mrad/s  or  m/s 

E,  is  the  spatial  frequency  in  cycles/mrad  or  cycles/m 

C.  OPTICAL  TRANSFER  FUNCTION  THEORY 

As  mentioned  in  the  previous  section,  the  optical  transfer  function  is  simply 
defined  as  the  Fourier  Transform  of  the  point  spread  function.  The  optical  transfer 
function  is  multiplied  with  the  input  spectrum  in  order  to  calculate  the  output  spectrum  of 
a  system  meaning  that  each  frequency  component  of  the  input  spectrum  is  weighted  and 
altered  by  a  factor  detennined  by  the  optical  transfer  function  [Ref.  9]. 

In  the  case  of  a  cascaded  system,  the  product  of  each  sub-system  transfer  function 
(assuming  the  sub-system  impulse  responses  being  linear  and  independent  of  each  other) 
would  detennine  the  overall  system  transfer  function.  Notice  that  it  is  more  practical  to 
do  the  computations  in  the  spatial  frequency  domain  than  to  do  them  in  the  spatial 
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domain  because  the  simple  product  of  the  transfer  functions  will  replace  the  complicated 
operation  of  the  multiple  convolution  integrals  in  the  spatial  domain  [Ref  10]. 

The  transfer  function  for  optical  systems,  in  other  words,  the  optical  transfer 
function  (OTF)  is  a  complex  quantity  as  mentioned  earlier  in  this  chapter.  The  magnitude 
of  OTF  yields  the  modulation  transfer  function  (MTF),  which  is  a  sine  wave  amplitude 
response  function.  By  convention,  MTF  approaches  zero  at  the  system  cut-off  spatial 
frequency,  which  indicates  the  highest  detail  of  the  target  that  can  be  reproduced  by  the 
system.  The  phase  of  the  OTF  gives  the  phase  transfer  function  (PTF)  that  denotes  the 
phase  shifts  brought  about  by  the  system  [Ref.  10].  Equation  (2.7)  obviously  shows  the 
relationship  between  OTF,  MTF  and  PTF  functions. 

OTF  ,  ij)  =  MTF  (£  ,  77  )e  JPTF  (2.7) 

The  effect  of  PTF  can  be  neglected  if  the  impulse  response  is  a  symmetric 
function.  In  that  case,  MTF  becomes  the  primary  physical  parameter  in  the  specification 
of  the  thermal  imaging  system  and  its  components.  Moreover,  in  the  areas  of  the  system 
requirement  analysis,  optimum  system  design  and  the  system  trade-off  analysis,  it  is  one 
of  the  system  parameters  most  commonly  manipulated  by  the  system  designers  and 
system  engineers. 

However,  using  MTF  theory  requires  some  assumptions.  Although  the  linear 
theory  and  the  MTF  theory  provide  very  practical  tools  to  the  system  evaluation  analysis, 
there  is  a  great  concern  about  the  assumptions  that  are  made  because  they  are  frequently 
violated  to  some  extent  in  real  system  applications.  For  example,  individual  detector 
responses  in  a  detector  array  being  uniform,  the  electronic  processing  being  linear  and  the 
system  mapping  being  single  valued  and  non-noisy  are  some  of  the  primary  violations 
that  have  been  experienced  [Ref.  10]. 

D.  RESOLUTION  AND  SENSITIVITY 

The  thennal  imaging  systems’  image  quality  specifications  and  sensor  trade-off 
studies  are  after  discussed  in  terms  of  spatial  resolution  and  sensitivity.  Resolution  refers 


11 


to  how  well  a  sensor  can  perceive  and  reproduce  the  smallest  detail  of  the  target  while 
sensitivity  gives  a  visual  measure  of  how  well  the  sensor  can  detect  a  target  signal  in  the 
presence  of  noise  [Ref.  4]. 

Although  resolution  and  sensitivity  are  conflicting  terms  most  of  the  time,  change 
in  one  affects  the  other.  The  relationship  between  sensitivity  and  resolution  is  indicated 
schematically  in  Figure  2.3,  which  reveals  that  there  is  a  significant  trade-off  between 
resolution  and  sensitivity.  Neither  resolution  nor  sensitivity  defines  the  overall  sensor 
performance  very  well,  but  together  they  produce  a  more  comprehensive  characterization 
[Ref.  9]. 
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Figure  2.3.  Resolution  and  Sensitivity  Relationship  [Ref.  9]. 

Sensitivity  becomes  especially  important  when  viewing  or  detecting  the  signals 
from  the  target  of  interest  immersed  in  the  noise  generated  by  a  variety  of  sources.  These 
sources,  either  environmental  or  sensor-based,  introduce  noise  in  every  stage  of  the 
thermal  imaging  process  [Ref.  9].  Therefore,  radiometry  which  describes  the  amount  of 
energy  leaving  the  object  space  and  arriving  at  the  sensor,  light  gathering  properties  of 
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optics  design,  detector  perfonnance  specifications  such  as  responsivity,  and  system  noise 
become  essential  tools  in  the  system  sensitivity  analysis  [Refs.  4  and  9]. 

Sensitivity  is  usually  described  by  a  noise  equivalent  parameter,  frequently 
defined  as  the  target-to-background  radiance  difference  giving  a  signal  to  noise  ratio 
(SNR)  of  one  at  the  sensor  (i.e.  NEI).  Noise  equivalent  temperature  difference  (NETD)  is 
a  commonly  used  noise  equivalent  parameter,  which  is  a  good  measure  of  sensitivity, 
measured  at  the  amplified  detector  output  of  the  thermal  imager.  [Refs. 9  and  10] 

On  the  other  hand,  resolution  is  the  capability  of  the  sensor  that  makes  the 
discrimination  tasks  of  detection,  recognition  and  identification  applicable,  when  a 
sufficient  level  of  sensor  sensitivity  is  achieved.  Resolution  depends  on  the  sensor  sub¬ 
system  parameters  such  as  optical  aberration  and  diffraction,  detector  size,  detector  shape 
and  focal  length  of  the  optical  aperture.  Moreover,  the  electronic  bandwidth,  the  display 
resolution  and  many  other  parameters  influence  the  sensor  overall  resolution  capacity. 
Notice  that  resolution  does  not  involve  the  noise  effects  [Refs.  4  and  9]. 

The  resolution  requirement  of  a  sensor  is  determined  by  the  most  stringent 
demand  of  the  recognition  task.  More  stringent  discrimination  tasks,  such  as 
identification,  are  mostly  performed  in  a  narrow  field-of-view.  Field-of-view  is  defined 
as  the  system’s  total  angular  space  in  the  scene  from  where  the  radiation  is  received.  Due 
to  the  technological  limitations  of  the  number  of  detectors  in  a  focal  plane  array  and  the 
number  of  corresponding  samples  taken  across  the  target,  a  narrow  field-of-view  provides 
more  samples  across  a  smaller  area  within  the  scene,  providing  more  details  on  the  target. 
Today  there  are  sensors  having  both  narrow  field-of-view  capability  providing  better 
resolution  and  wide  field-of-view  providing  better  coverage  over  large  areas  [Ref.  9]. 

E.  SAMPLING  THEORY  AND  ALIASING 

Sampling  is  an  unavoidable  process  in  all  real  thermal  imaging  systems  in  which 
an  array  of  the  samples  is  taken  from  a  two-dimensional  image  function  at  a  discrete  set 
of  points  in  the  image  plane.  These  samples  will  represent  the  original  image  function  if 
the  samples  are  taken  sufficiently  close  to  each  other.  Consequently,  the  original  function 
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would  be  reconstructed  by  a  simple  interpolation  applied  to  the  sampled  values  with 
considerable  accuracy  to  form  the  original  image  [Ref.  12]. 


Horizontal  Pixel  Pitch 


Vertical  Pixel  Pitch 


Figure  2.4.  Spatial  Sampling  in  a  FPA. 

Sampling  can  be  performed  in  several  ways  and  in  several  dimensions.  In 
scanning  systems,  the  detector  acts  as  a  sampling  aperture  while  the  analog-to-digital 
(A/D)  converter  acts  as  a  second  sampler  in  case  the  detector  output  is  digitized.  In 
staring  systems,  the  discrete  detectors  present  in  a  focal  plane  array  sample  the  scene 
spatially  [Ref.  4].  One  primary  characteristic  of  the  staring  systems  is  that  they  introduce 
sampling  in  two  dimensions.  The  spatial  sampling  process  in  a  focal  plane  array  is  shown 
above  in  Figure  2.4.  The  outer  rectangular  frame  represents  the  FPA  in  which  the  discrete 
detector  elements  are  located.  The  detector  elements  are  assumed  to  be  rectangular 
because  it  is  more  convenient  for  analysis  and  it  is  the  most  common  geometry 
encountered.  The  spatial  sampling  frequency  depends  on  the  distance  between  two 
adjacent  detector-element  centers,  called  the  pixel  “pitch”.  Therefore,  spatial  sampling 
period  is  calculated  as: 


T 
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where  T  is  the  sampling  period  in  radians,  Ad  is  the  pixel  pitch  in  meters  and  /optics  is  the 
focal  length  of  the  sensor  optics  in  meters.  The  sampling  frequency  is  defined  as: 


fs 


m 


(2.9) 


where  fs  is  the  sampling  frequency  in  samples  (cycles)  per  radian.  It  is  clear  that  the 
sampling  frequency  in  the  vertical  direction  is  equal  to  the  sampling  frequency  in  the 
horizontal  direction  for  square  detector  geometry. 


The  sampling  theory  requires  that  if  the  original  image  being  sampled  contains 
information  with  frequencies  above  half  the  sampling  frequency,  which  is  called  the 
Nyquist  frequency,  the  final  image  will  contain  spurious  information  after  the 
reconstruction  process.  In  this  process,  a  frequency  component,  which  is  above  the 
Nyquist  limit  in  the  original  image,  appears  as  a  lower  frequency  component  after 
sampling  due  to  the  phenomenon  called  “aliasing”.  The  Nyquist  frequency  (f/)  is  defined 
as: 


L 


(A) 
U ) 


(2.10) 


The  result  of  aliasing  can  be  observed  as  distortions  and  degradation  in  the  image 
produced. 

One-dimensional  frequency  spectra  of  an  original  signal  before  and  after  sampling 
at  differing  sampling  rates  are  given  in  Figure  2.5.  The  band-limited  signal  in  Figure 
2.5(a),  which  is  also  called  the  baseband,  is  the  original  signal  before  the  sampling 
process.  It  is  important  to  notice  that  sampling  theory  can  only  be  applied  to  band-limited 
functions  for  an  exact  reconstruction  process  of  the  original  functions. 

The  replication  property  of  delta  function  suggests  that  sampling  of  a  signal  in  the 
spatial  domain  corresponds  to  the  replication  of  the  band-limited  original  spectrum  in  the 
spatial  frequency-domain.  As  a  result  of  this  concept,  the  sidebands  represent  the 
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replications  of  the  original  signal  at  frequencies  Kfs  after  being  sampled,  where  K  is  an 
integer  (K  =  1,2,3...).  If  the  sampling  frequency  were  increased  so  that  there  would  be  no 
overlapping,  the  sampled  image  would  be  perfectly  recovered.  This  phenomenon  is  called 
oversampling  and  is  depicted  in  Figure  2.5(b).  Oversampling  does  not  necessarily  mean 
that  there  is  excessive  sampling  in  the  process. 

Figure  2.5(c)  shows  the  critical  sampling  case  in  which  the  sampling  frequency  is 
arranged  so  that  the  Nyquist  frequency  requirement  is  barely  reached  (fs  =2fi1).  Only  a 
perfect  low-pass  filter  can  recover  the  original  input  image  by  removing  the  sidebands. 

If  a  frequency  of  the  original  signal  is  greater  than  the  Nyquist  frequency  fy  which 
is  the  half  of  the  sampling  rate,  the  system  is  undersampled.  The  overlap  of  the  baseband 
and  the  sidebands  causes  aliasing  as  in  Figure  2.6  where  /[,  =/s  =2  /N.  There  would  be  an 
ambiguity  within  the  overlapping  frequency  band  because  it  would  be  impossible  to 
know  the  origin  of  the  reconstructed  spectrum  frequency  components  [Refs.  4  and  11]. 
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Figure  2.5.  Sampling  Process  at  Different  Sampling  Frequencies  (A)  Original 
Band-Limited  Signal  (B)  Frequency  Spectrum  after  Sampling  and  (C)  when 
(fs  =2 /h),  the  Bands  Adjoin  [Ref.  4]. 
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RELATIVE  FREQUENCY 


Figure  2.6.  The  Overlap  of  the  Baseband  and  the  Sidebands  Causes  Distortion  in 
the  Output  Spectrum  due  to  Undersampling  [Ref.  4], 

Many  methods  have  been  suggested  to  reduce  the  aliasing  effect  in  the 
final  image.  One  of  them  suggests  equating  the  pitch  to  the  sampling  period  using  a  fill 
factor  of  100%  where  fill  factor  is  defined  as  the  ratio  of  the  FPA  active  area  to  the  FPA 
total  area.  Flow  ever,  100%  fill  factor  in  FPAs  is  inapplicable  due  to  the  fabrication 
technology  limitations  today.  Another  method  suggests  increasing  the  detector 
dimensions  decreasing  the  spectral  width,  but  that  reduces  the  modulation  of  detectable 
spatial  frequencies.  Consequently,  there  seems  to  be  a  trade-off  between  system  MTF  and 
the  aliasing  phenomenon.  An  alternate  way  might  be  to  increase  the  sampling  rate.  This 
can  be  accomplished  either  by  increasing  the  number  of  detector  elements  in  an  array  or 
by  decreasing  the  field-of-view.  Increasing  the  number  of  detector  elements  is  difficult 
again  due  to  the  technological  limitations,  while  decreasing  the  field-of-view  approach 
might  be  undesirable. 

The  last  but  the  most  popular  method  being  suggested,  which  is  also  in  use 
today,  is  the  micro-scanning  process.  This  technique  is  a  way  of  increasing  the  sampling 
rate  by  moving  the  image  around  on  the  detector  array  by  fractions  of  a  pixel  size  and 
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recording  the  successive  images  for  each  shift.  The  recorded  images  are  then  combined  to 
produce  a  final  image.  Although  the  micro-scanning  technique  is  applicable  to  today’s 
thermal  imagers,  it  brings  the  trade-off  between  aliasing  or  integration  time  together  [Ref. 
13]. 

F.  SYSTEM  PERFORMANCE  PARAMETERS 

The  following  section  describes  the  three  basic  parameters:  Noise  Equivalent 
Temperature  Difference  (NETD),  Minimum  Resolvable  Temperature  Difference 
(MRTD)  and  Minimum  Detectable  Temperature  Difference  (MDTD)  that  are  used  in  the 
performance  specifications  of  the  thermal  imaging  systems. 

1.  Noise  Equivalent  Temperature  Difference  (NETD) 

NETD  is  an  intennediate  sensitivity  parameter  defined  as  the  target-to- 
background  differential  temperature  generating  a  peak  signal  to  rms  noise  ratio  of  one 
while  imaging  relatively  large  targets  [Ref.  14],  For  system  measurements,  NETD  is 
determined  at  the  analog  video  or  at  the  output  of  the  monitor  through  a  standard  filter, 
which  corresponds  to  the  filtering  performed  by  the  electronics  prior  to  the  measurement 
point  [Refs.  15  and  16].  Basically,  NETD  for  a  typical  thermal  imaging  system  is 
determined  by  the  simple  equation  given  below: 


NETD  = 


'  AT  ' 


(2.11) 


where  AT  represents  the  temperature  difference  between  the  target  and  the  background. 
Notice  that  when  the  denominator  value  is  set  to  one,  NETD  becomes  equal  to  AT  by 
definition. 


For  a  scanning  system,  Shumaker  [Ref.  14]  defines  NETD  with  the  equation 
given  below: 


NETD  = 


\Q{FOVxFOVyFrNosNJ12 
(tzN Drjsc)112  DAxAyD**rjcs  (dN  /  dT)r0 


(2.12) 
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where: 


Nos  is  the  over-scan  ratio 

Nss  is  the  serial  scan  ratio 

Nd  is  the  number  of  detectors 

8N/8T  is  the  thennal  gradient  (“the  derivative  of  Planck’s  equation”  in  watt  cm'2 
K'1  sr'1) 

T)sc  is  the  scan  efficiency 

Tjcs  is  the  cold-shield  efficiency 

D  is  the  optics  diameter  (m) 

D  is  the  band  averaged  detector  detectivity  (cm  Hz'  W'1) 

Ax  is  the  in-scan  detector  angular  subtense  (mRad) 

Ay  is  the  cross-scan  detector  angular  subtense  (mRad) 

FOVx  is  the  in-scan  lield-of-vicw  (mRad) 

FOVy  is  the  cross-scan  field-of-view  (mRad) 

Fr  is  the  frame  rate  (1 /sec) 

To  is  the  optical  transmittance 

With  the  growth  of  the  thennal  imaging  technology,  NETD  has  become  a  major 
performance  predictor  parameter.  Its  use  simplified  the  formulations  for  the  more 
comprehensive  prediction  parameters  such  as  MRTD  and  MDTD.  However,  NETD  does 
not  include  the  temporal  and  spatial  integration  effects  of  the  eye,  and  characterizes  only 
the  temporal  detector  noise,  instead  of  accounting  for  overall  system  noise  from  various 
sources. 

2.  Minimum  Resolvable  Temperature  Difference  (MRT) 

Minimum  resolvable  temperature  difference  (MRTD  or  MRT)  is  the  primary 
summary  measure  of  a  thermal  imaging  system  performance  including  resolution, 
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sensitivity,  and  the  observer  [Ref.  17].  It  is  defined  as  the  temperature  difference  between 
a  uniform  background  and  the  bars  of  a  four-bar  target,  which  makes  the  individual  bars 
just  resolvable  through  a  thennal  imager  to  a  trained  observer  [Ref.  11].  Four-bar  patterns 
are  standardized  by  a  7: 1  aspect  ratio,  which  makes  the  overall  pattern  a  square. 

A  typical  MRT  curve  can  be  formed  with  respect  to  the  relative  sizes  of  the  bar 
patterns  as  in  Figure  2.7.  Some  of  the  substantial  characteristics  of  MRT  can  be  observed 
in  this  figure.  First,  MRT  is  not  a  single  number;  furthermore,  it  is  a  function  of  spatial 
frequency  [Ref.  14],  Second,  the  temperature  difference  to  resolve  the  four-bar  chart 
increases  as  the  bar  dimensions  decrease  [Ref.  11].  The  MRT  curve  ultimately  shows  an 
asymptotic  characteristic  at  a  spatial  frequency  around  the  1/DAS  value,  which  makes  the 
TIS  to  be  limited  by  the  detector  characteristics,  since  even  very  large  temperature 
differences  will  not  make  the  individual  bars  resolvable  [Ref.  14].  DAS  is  the  detector 
angular  subtense,  which  describes  the  achievable  sensor  resolution  relative  to  the  detector 
dimension  limitations.  Dividing  the  detector  dimension  (width  or  height)  by  the  optical 
focal  length  yields  DAS  in  that  dimension  [Ref.  9].  A  final  aspect  of  MRT  to  be  noticed 
is  that  it  is  not  an  objective  measurement  because  it  comprises  the  human  eye-brain 
performance. 
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Figure  2.7.  A  Typical  MRT  Plot  [After  Ref.  11]. 


Numerous  equations  with  various  modifications  and  approximations  have  been 
proposed  for  MRT  prediction.  One  general  form  for  a  scanning  system  is  given  by 
Shumaker  [Ref  14]  in  the  following  formula: 


MRT(v )  = 


2 0 SNR T [fO Vx FO Vy v 2 )'  2  plJ2 
t0DD** {xNDrjsc  )1/2  ncs  (AxAv)1/2M7F5  (v\Lte  f2  dN/dT 


(2.13) 


Where: 

v  is  the  spatial  frequency  (cycles/mRad) 
te  is  the  eye  integration  time  (s) 
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8N/8T  is  the  thermal  gradient  (“the  derivative  of  Planck’s  equation”  in  watt  cm'2 
K'1  sr'1) 

Nd  is  the  number  of  detectors 

r)sc  is  the  scan  efficiency 

r)cs  is  the  cold-shield  efficiency 

FOVx  is  the  in-scan  field-of-view  (mRad) 

FOVy  is  the  cross-scan  field-of-view  (mRad) 

Ax  is  the  in-scan  detector  angular  subtense  (mRad) 

Ay  is  the  cross-scan  detector  angular  subtense  (mRad) 

MTFs  is  the  overall  system  modulation  transfer  function 
SNRT  is  the  threshold  signal-to-noise  ratio 
D  is  the  diameter  of  the  optical  aperture  (m) 

D**  is  the  band-average  detector  specific  detectivity  for  271  steradians  (cm  Hz/2 
W'1) 

r0  is  the  optical  transmittance 

L  is  the  length-to-width  aspect  ratio  of  the  target  bar 
px  is  the  noise  filter  factor 

3.  Minimum  Detectable  Temperature  Difference  (MDT) 

Although  the  definition  of  the  minimum  detectable  temperature  difference 
(MDTD  or  MDT)  is  almost  the  same  as  the  definition  of  the  MRT,  the  difference 
between  the  two  arises  from  the  representation  of  the  target.  MDT  is  a  function  of  target 
size  where  the  target  is  represented  by  an  isolated  square.  Therefore,  MDT  is  defined  as 


23 


the  differential  temperature  between  the  square  target  and  a  uniform  background,  which 
makes  the  square  just  detectable  through  a  thermal  imager  to  a  trained  observer  [Ref  14]. 

A  representative  MDT  curve  can  be  formed  with  respect  to  the  relative  sizes  of 
the  square  targets  as  in  Figure  2.8.  Similarly,  MDT  is  plotted  as  a  function  of  spatial 
frequency,  which  can  be  calculated  as  the  reciprocal  of  twice  the  angular  subtense  of  the 
side  of  the  square  target  [Ref.  14].  The  plot  shows  that  the  required  temperature 
difference  to  detect  the  square  target  increases  as  the  target  size  decreases  and  the  MDT 
curve  does  not  show  an  asymptotic  characteristic  as  the  MRT  curve  does  because  even  a 
point  target  can  be  made  detectable  by  increasing  its  temperature. 
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Figure  2.8.  A  Typical  MDT  Plot  [Ref.  11]. 

A  commonly  known  equation  for  MDT  is  given  by  Shumaker  [Ref.  14]  as: 
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(2.14) 


MDT(v) 


20SNRT(FOVxFOVvAxAy)12  (qt  +  r2 ) 
7TT0DAxAyD**  8N/8T  rjcs  ( NDrjJe  )1/2QT  (rs2  +  r2B  +  Q  T  f2 


where: 

v  is  the  spatial  frequency  (cycles/mRad) 
te  is  the  eye  integration  time  (s) 

8N/8T  is  the  thennal  gradient  (“the  derivative  of  Planck’s  equation”  in  watt  cm'2 
K'1  sr'1) 

Nd  is  the  number  of  detectors 

rise  is  the  scan  efficiency 

r|cs  is  the  cold-shield  efficiency 

FOVx  is  the  in-scan  field-of-view  (mRad) 

FOVy  is  the  cross-scan  field-of-view  (mRad) 

Ax  is  the  in-scan  detector  angular  subtense  (mRad) 

Ay  is  the  cross-scan  detector  angular  subtense  (mRad) 

rs  is  the  overall  resolution  of  the  imaging  system  (mRad) 

Hi  is  the  overall  resolution  of  the  back-end  components  including 

components  from  the  detector  electronics  through  the  observer(mRad) 

Qt  is  the  solid  angular  subtense  of  the  square  target  (mRad) 

SNRT  is  the  threshold  signal-to-noise  ratio 

D  is  the  diameter  of  the  optical  aperture  (m) 
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D**  is  the  band-average  detector  specific  detectivity  for  2 71  steradians  (cm  Hz  2 

W1) 

t0  is  the  optical  transmittance 

The  equation  above  indicates  a  weaker  dependence  on  the  system  MTF  when 
compared  to  the  MRT  prediction  characteristics.  On  the  other  hand,  it  is  directly  related 
to  the  detection  capability  of  a  noise-limited  system,  which  makes  MDT  a  measure  of 
thennal  system  sensitivity.  It  differs  from  NETD  in  that  the  observer  eye-brain 
performance  is  accounted  for  in  MDT  calculations  [Ref.  11]. 
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III.  FLIR-92  TIS  PERFORMANCE  MODEL 


The  purpose  of  this  chapter  is  to  summarize  the  basic  principles  of  the  current 
standard  FLIR-92  model  that  are  relevant  to  the  content  of  the  proceeding  chapters. 

FLIR-92  is  a  computer-based  model  used  as  a  system  evaluation  tool,  which 
predicts  standard  summary  performance  measures  for  thermal  imaging  systems.  It 
operates  in  the  DOS  environment  and  calculates  the  modulation  transfer  function  (MTF), 
noise  equivalent  temperature  difference  (NETD),  minimum  resolvable  temperature 
difference  (MRTD)  and  minimum  detectable  temperature  difference  (MDTD)  utilizing 
the  primary  sensor  parameters.  It  has  a  principal  function  of  predicting  whether  or  not  a 
system  achieves  the  desired  MTF,  system  noise,  MRTD  and  MDTD  detennined 
necessary  to  perfonn  mission-specific  target  discrimination  and  acquisition  tasks  [Ref.  2]. 

FLIR-92  can  be  used  in  modeling  both  scanning  and  staring  thermal  imagers  that 
operate  in  mid-infrared  (3-5  pm.)  and  far-infrared  (8-12  pm.)  spectral  bands.  However,  it 
cannot  be  used  for  any  type  of  electro-optical  sensors  other  than  thermal  imagers  (e.g. 
night  vision  devices).  Furthermore,  acquisition  and  discrimination  range  performance 
prediction  is  beyond  the  capabilities  of  FLIR-92  [Ref.  2]. 

The  details  of  the  model  will  be  covered  under  three  topics: 

•  Modulation  Transfer  Functions  (MTFs) 

•  3-D  Noise  Concept 

•  MRTD  and  MDTD  Predictions 

A.  MODULATION  TRANSFER  FUNCTIONS  (MTFS) 

The  FLIR-92  model  assumes  the  thennal  imaging  system  under  evaluation  to  be 
linear  to  enable  the  use  of  linear  system  analysis  in  calculations.  The  evaluation  process 
starts  with  recording  the  system  parameters  in  the  list  of  inputs,  which  will  be  used  in  the 
component  MTF  calculations.  At  this  point  the  analyst  has  the  flexibility  of  either  using 
the  built-in  default  system  parameters  or  configuring  them  as  desired  [Ref.  2], 


27 


In  FLIR-92,  system  component  MTFs  are  grouped  into  three  categories  for 
modeling  purposes:  pre-filters,  temporal  post-filters  and  spatial  post-filters.  The  MTF 
groups  and  their  components  are  listed  in  Table  3.1  and  the  corresponding  equations  used 
in  the  calculations  are  presented  in  Appendix  A.  As  a  result  of  the  linear  system 
approach,  FLIR-92  calculates  the  overall  system  transfer  function  by  simply  multiplying 
together  the  entire  set  of  component  MTFs. 


Pre-filter 

MTFs 

OPTICS  MTFS 

-Diffraction-limited  MTF 
-Geometric  Blur  MTF 

Detector  Spatial  MTF 

Focal  Plane  Array  Integration  MTF 

Sample-scene  Phase  MTF 

Image  Motion  MTFs 

-Linear  Image  Motion  MTF 
-Random  Image  Motion  MTF 
-Sinusoidal  Image  Motion  MTF 

Temporal  Post-filter 
MTFs 

Detector  Temporal  MTF 

Electronics  Low  Frequency  Response 

Electronics  High  Frequency  Response 

Boosting  MTF 

Spatial  Post-filter 
MTFs 

Electro-optical  Multiplexor  MTF 

Digital  Filter  MTF 

Display  MTF 

-CRT  Display  MTF 

CCD  Charge  Transfer  Efficiency  MTF 

Display  Sample  and  Hold  MTF 

Eye  MTF 

-Non-limiting  Eye  MTF 
-Limiting  Eye  MTF 

Table  3.1.  FLIR-92’s  Sub-system  MTFs  Grouped  into  Three  Categories. 

FLIR-92  gives  the  MTF  results  in  two  different  forms.  The  short  output  lists  the 
MTF  results  for  only  three  main  groups  (pre-filters  MTF,  temporal  post-filters  MTF  and 
spatial  post-filters  MTF),  while  the  long  output  gives  each  single  MTF  calculation 
selected  by  the  analyst  under  these  groups. 
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B.  3-D  NOISE  CONCEPT 

In  the  early  models,  the  system  noise  was  adequately  quantified  by  a  one¬ 
dimensional  NETD  value,  which  was  assumed  to  be  equal  to  the  predominant  detector 
noise  for  the  first-generation  thennal  imagers.  With  the  introduction  of  the  advanced 
scanning  and  staring  thennal  imagers,  temporal  detector  noise  can  no  longer  sufficiently 
characterize  the  system  noise  due  to  the  use  of  focal  plane  array  technology  and  the 
associated  non-uniformity  conection  schemes,  which  are  considered  as  complex  noise 
sources  for  the  sensors.  Moreover,  Time  Delay  and  Integration  (TDI),  digital  sampling 
and  digital  processing  techniques,  which  were  used  to  increase  the  signal-to-noise  ratio 
(SNR),  introduced  additional  noise  terms  to  the  system  [Ref.  18]. 

It  was  observed  that  the  noise  patterns  produced  by  these  thermal  imagers  show 
direction  dependent  characteristics.  In  order  to  provide  more  accurate  predictions, 
D’Agostino  and  Webb  developed  a  new  approach  that  takes  the  new  noise  components 
into  consideration.  Their  methodology  suggests  dividing  the  total  noise  present  at  the 
output  of  the  sensor  into  a  group  of  seven  direction  dependent  noise  components,  which 
simplifies  the  analysis,  the  understanding  and  also  the  incorporation  of  the  noise  terms 
into  model  formulations.  Each  noise  component  shows  specific  characteristics  in  each 
direction  in  a  three-dimensional  coordinate  system  that  is  shown  in  Figure  3.1  [Ref.  2]. 


Figure  3.1,  Three-dimensional  Coordinate  System  [Ref.  2], 
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The  script  t  represents  the  temporal  (time)  dimension  or  sequential  frames,  while 
h  represents  the  horizontal  direction  or  the  sequential  columns  within  a  frame  and  v 
represents  the  vertical  direction  or  the  sequential  rows  within  a  frame.  The  horizontal 
and  vertical  directions  provide  the  spatial  information  in  the  original  image.  Depending 
on  the  type  of  the  thennal  imager,  the  horizontal  direction  may  represent  a  spatial 
location  in  a  staring  sensor  or  time  in  a  scanning  system  [Ref.  4]. 


Noise 

Description 

Source 

OXVH 

Random  spatio-temporal  noise 

Basic  detector  temporal  noise 

<3XV 

Temporal  row  noise  (i.e.  line 
bounce) 

Line  processing, 

1/f,  readout 

<3th 

Temporal  column  noise  (i.e. 
column  bounce) 

Scan  effects 

C?VH 

Random  spatial  noise  (i.e.  bi¬ 
directional  fixed  pattern  noise 

Pixel  processing,  detector-to- 
detector  non-uniformity,  1/f 

av 

Fixed  row  noise  (i.e.  line-to-line 
non-uniformity) 

Detector-to-detector  non¬ 
uniformity 

oh 

Fixed  column  noise  (i.e.  column  - 
to-  column  non-uniformity) 

Scan  effects,  detector-to-detector 
non-uniformity 

CTt 

Frame-to-frame  noise  (i.e.  frame 
bounce) 

Frame  processing 

S 

Mean  of  all  noise  components 

Table  3.2.  3-D  Noise  Component  Descriptions  [Ref.  2], 

The  components  that  define  the  overall  noise  are  measured  at  the  system  output 
and  they  are  finally  converted  to  temperature  (degrees  C)  in  the  same  manner  as  done  for 
NETD.  The  components,  their  descriptions  and  sources  are  as  given  in  Table  3.2.  The 
subscripts  indicate  the  directions  in  which  the  noise  fluctuations  about  a  mean  of  zero 
occur.  The  missing  subscripts  indicate  the  directions  in  which  the  averaging  operation  is 
done.  The  averaging  operation,  which  is  done  using  the  “D  operator”  within  the 
coordinate  system,  enables  the  decomposition  of  the  complex  noise  into  its  components. 
A  detailed  explanation  of  the  D  operator  and  its  function  is  provided  in  Appendix  B. 
Averaging  complex  noise  in  a  specific  direction  yields  the  desired  noise  component 
present  in  the  other  directions  in  which  no  averaging  operation  is  done.  For  example,  axv 
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represents  the  rms  noise  value  that  is  present  in  the  temporal  and  vertical  directions 
calculated  by  averaging  the  complex  noise  in  the  horizontal  direction  [Refs.  2  and  9]. 

The  noise  components  are  assumed  to  be  independent  of  each  other.  As  a  result  of 
this  the  total  system  noise  is  defined  as  the  root  sum  square  of  the  noise  components.  This 
relation  is  shown  in  the  following  equation: 


Q  =  (<Jtvh  +  <J 


TV 


2 

CJ  TH 


2  2 

CT  VH  ^  CJ  V 


2 

CJh 


A0-5 

+  OV 


(3.1) 


where  Q  represents  the  overall  system  noise. 

All  noise  components  may  not  be  seen  in  every  imaging  system  and  depending  on 
the  specific  system  design  one  or  more  may  dominate  the  others.  They  can  be  determined 
either  from  real  system  measurements  or  from  estimates  except  for  the  random  spatio- 
temporal  noise  component  (cttvh).  otvh  is  the  basic  detector  noise  that  is  seen  due  to 
arriving  photons  (shot  noise)  and  flowing  current  (Johnson  noise).  Therefore,  it  can  be 
calculated  in  a  similar  fashion  to  that  used  in  the  first  generation  systems  and  can  be 
related  to  the  NETD  value  as  seen  in  the  following  equation. 


C 7 


TVH 


=  NETD  x 


(3.2) 


where  A/\-  is  the  equivalent  noise  bandwidth  for  NETD  calculation  and  Afp  is  the  actual 
system  noise  bandwidth  associated  with  the  electronics  prior  to  the  display  [Ref.  2].  It  is 
obvious  from  the  equation  above  that  FLIR-92  calculates  cjtvh  using  the  actual  system 
bandwidth  rather  than  using  an  artificial  reference  bandwidth  as  done  for  NETD. 

The  rest  of  the  noise  components  other  than  cttvh  can  be  classified  into  two 
groups.  Random  spatial  noise  or  so-called  fixed  pattern  noise  (<tvh),  fixed  row  noise  (av) 
and  fixed  column  noise  (cth)  are  all  time  independent  noise  components  [Ref.  19].  They 
occur  due  to  the  detector  and  electronics  non-uniform  responses  over  the  focal  plane 
array  that  might  originate  from  the  detector  material  variations,  impurities  or 
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imperfections,  electrical  circuit  tolerances,  and  so  on  [Ref.  9].  The  time  dependent  noise 
components,  which  can  be  listed  as  frame-to-frame  noise  (aT),  temporal  row  noise  (arv) 
and  temporal  column  noise  (gth)  are  mostly  seen  in  scanning  imaging  systems,  ax  is 
assumed  negligible  so  that  it  is  not  included  in  the  model  predictions  [Ref.  19].  The 
global  average,  the  mean  of  all  noise,  is  denoted  by  the  symbol  S. 

As  mentioned  previously  the  observer  is  one  of  the  system  components  in  FLIR- 
92  predictions.  To  incorporate  the  observer’s  eye-brain  effects,  FLIR-92  uses  a 
“synchronous  integrator  model”  for  MRTD  predictions  while  it  uses  the  “matched  filter 
model”  for  the  MDTD  predictions.  The  eye-brain  effects  are  included  into  the  three- 
dimensional  noise  methodology  using  the  eye-brain  integration  factors  ET,  EH,  Ev,  which 
describe  the  temporal  integration,  horizontal  integration  and  vertical  integration 
respectively  [Ref.  4].  The  integration  factors  are  effective  only  when  the  integration 
direction  and  the  noise  component  direction  coincide.  As  a  consequence  of  this  approach 
the  total  system  noise  term  can  be  modified  as  the  root  sum  square  of  the  components  as 
in  Equation  3.3. 

U  (fs)  =  (<t2tv„EtEVz  (fs)E„t  (  fs )  +  a\«EVz  (  fs  )E„z  (  fs )  +  a2T„ETEHz  (  fs )  + 
<J2t,EtEVz  (fs)  +  a2rETE,,z  (fs)  +  <72hEHz  (/s)A  «-3' 

Note  that  the  gt  component  has  been  neglected  in  the  equation  above  [Ref.  2]. 

C.  MRTD  AND  MDTD  PREDICTIONS 

FLIR-92  uses  different  eye-brain  models  for  the  predictions  of  MRTD  and 
MDTD  as  mentioned  in  the  previous  section.  For  predicting  MRTD  it  uses  the 
synchronous  integrator  model,  which  assumes  the  eye-brain  components  to  be  a  spatially 
integrating  system  over  the  image,  while  ignoring  blurring  effects  caused  by  the  finite 
sized  sensor  apertures.  On  the  other  hand,  for  the  prediction  of  MDTD  it  uses  the 
matched  filter  model,  which  relates  the  eye-brain  effects  to  a  matched  filter  resulting  in 
the  increase  of  SNR  ratio.  Both  models  yield  the  same  results  for  the  periodic  targets  such 
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as  four-bar  pattern.  Due  to  the  simplicity  of  integrating  the  algorithm  into  FLIR-92,  the 
synchronous  integrator  model  is  preferred  for  the  MRTD  calculations  [Ref.  2]. 

Once  all  noise  components  are  determined  either  by  measurement  or  calculation, 
the  components  in  the  horizontal  direction  are  integrated  into  the  horizontal  MRTD 
equation  while  the  components  in  the  vertical  direction  are  integrated  into  the  vertical 
MRTD  equation  in  the  same  way.  Including  the  noise  effects  degrades  the  overall  system 
performance.  A  general  interpretation  of  the  MRTD,  which  clearly  shows  the  inclusion  of 
the  correction  factor  into  the  MRTD  prediction,  is  presented  in  the  equation  below. 


MRTDz  (fs) 


ff  K1 


SNRth  crTVH  kz  (fs ) 
MTFSYS(fs ) 


ETEVz(fs)EHz(fs)] 


12 


(3.4) 


The  subscript  z  in  the  equation  refers  to  the  direction  of  interest  either  horizontal 
or  vertical.  The  other  parameters  and  their  descriptions  are  given  in  Table  3.3. 


PARAMETERS 

DESCRIPTION 

SNRthr 

Threshold  signal-to-noise  ratio 

cttvh 

Random  spatio-temporal  noise 

kz 

Noise  correction  function 

MTFsys 

Overall  system  modulation  transfer  function 

Et 

Eye  temporal  integration  factor 

EHz(fs) 

Eye  horizontal  integration  factor  in  the  direction  of  interest 

Evz(fs) 

Eye  vertical  integration  factor  in  the  direction  of  interest 

Table  3.3.  MRTD  Equation  Parameters  [Ref.  2], 

In  the  equation  above,  aTvH  (the  random  detector  noise  or  spatio-temporal  noise) 
is  the  only  term  that  has  units  (degrees  Celsius).  The  rest  of  the  terms,  which  are 
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dimensionless,  are  functions  of  spatial  frequency.  The  U.S.  Army  Night  Vision  and 
Electronic  Sensors  Directorate  (NVESD)  recommends  a  value  of  2.5  for  the  threshold 
SNR  and  a  value  of  0.1  seconds  for  the  eye  integration  time  Et  as  reasonable  averages  for 
those  quantities  [Ref.  2]. 


The  gtvh  noise  tenn,  which  is  sometimes  approximated  with  NETD,  is  generally 
given  by: 


CJr 


f  X2 


sw 


7tt0^\  \D\A,300)— (A)dA 


300 


where: 

fi NO 

A  fP 


'to 

Ad 

D*(k,300) 

(8W/8T) 


is  the  optical  f-number 
is  the  system  noise  bandwidth 
is  the  optical  transmittance 
is  the  detector  area 

is  the  detector  noise-limited  spectral  detectivity 
is  the  thermal  derivative  of  Planck’s  Law  [Ref.  2] 


(3.5) 


The  noise  correction  function  for  the  horizontal  direction  is  given  as: 
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K(fs)  = 


1  + 


<7  vh  i  (7  m 
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<72h 


.  1/2 


1 


<72tvh  Er  <J~tvh  E v  (fs)  CT2tvh  ETE  (fs) 


(3.6) 


The  correction  function  for  the  vertical  direction  is  given  as: 


K(fs)  = 


i  + 


^2  i  ^2 

(T  VH  1  <7  TV 


1 


<72v 


.  1/2 


1 


<72tvh  Et  <7~tvh  Eh  (fs)  <72tvh  ETEHy(fs)j 


(3.7) 
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It  is  important  to  note  that  the  eye-brain  integration  factors  and  the  noise 
components  in  one  direction  can  only  modify  the  correction  function  in  that  direction. 
The  eye-brain  integration  factors  are  calculated  differently  for  staring  and  scanning 
systems  except  for  the  time  integration  factor,  which  is  given  in  the  equation  below. 

E=-°^-  (3.8) 

Frte 

where: 

cct  is  the  temporal  correlation  factor 

Fr  is  the  system  frame  rate  (Hz) 

Te  is  the  eye  integration  time  (s). 

The  spatial  integration  factors  for  staring  and  scanning  systems  are  calculated  as 
in  the  equations  given  in  Appendix  C. 


spatial  frequency 

Figure  3.2.  A  Typical  2-D  MRTD  Plot. 
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The  FLIR-92  model  calculates  the  horizontal  and  the  vertical  MRTD  values  for  a 
subset  of  spatial  frequencies.  Furthermore,  the  model  calculates  the  two-dimensional 
MRTD  values  and  their  corresponding  spatial  frequencies  from  the  geometric  mean  of 
the  horizontal  and  vertical  MRTD  spatial  frequency  values  as  presented  in  Figure  3.2. 
Since  each  MRTD  component  is  weighted  with  respect  to  its  frequency  axis,  the  MRTD2d 
curve  asymptotes  to  the  geometric  mean  value  of  the  horizontal  and  vertical  spatial  cutoff 
frequencies  as  seen  in  the  figure  below. 

The  FLIR-92  model  predicts  MDTD  using  Equation  3.9,  which  is  very  similar  to 
the  MRTD  equation. 


MDTD  (fs) 


SNR  m  cr  TVH  k  MDT  (fs ) 
l  ATQH(fs)Qv(fs)  J 


(ETEH(fs)Ev(fs)) 


1/2 


(3.9) 


fs  is  the  spatial  frequency  in  units  of  mrad'1  where  AT  is  the  target  area.  The  overall 
system  MTFs  in  the  horizontal  direction  and  in  the  vertical  direction  are  embedded  in  the 
quantities  Qn(fs)and  Qv(fs)respectively. 


QH(fs)  is  calculated  as: 
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Qv(fs)  is  calculated  as: 


QAfs)=\H;rs,(v) 


sin 


7IV- 


\  JsJ 


fs 


7TV- 


fs 


dv 


(3.11) 


36 


In  the  two  equations  given  above,  Hsys  represents  the  overall  system  MTF  in  the 
given  direction  and  the  sine  tenn  ((sinX)/X)  in  the  brackets  represents  the  Fourier 
Transform  of  the  square  target  by  definition.  For  both  staring  and  scanning  systems,  the 
same  equations  are  used  in  MDTD  predictions  [Ref.  2].  The  other  terms  in  the  MDTD 
formula  are  calculated  using  the  same  equations  as  defined  in  the  MRTD  case  except  for 
the  correction  function  kMDi{fs)  given  in  Equation  3.12. 
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IV.  THE  VISIBILITY  MODEL 


The  main  purpose  of  this  chapter  is  to  overview  the  basic  principles  of  the 
Visibility  concept,  which  was  proposed  by  R.J.Pieper  and  A.W.Cooper  of  the  Naval 
Postgraduate  School  in  1994.  This  chapter  will  also  serve  as  an  introduction  to  the  Virtual 
Thermal  Image-Processing  model  that  will  be  presented  in  Chapter  VI  because  the 
Virtual  Thermal  Image-Processing  model  is  based  on  the  Visibility  model.  A  short 
review  of  the  latest  modifications  in  the  Visibility  model  is  included  here  because  these 
amendments  are  also  utilized  in  constructing  the  Virtual  Thermal  Image-Processing 
model.  The  chapter  will  close  with  the  presentation  of  the  objective  MRTD  predicting 
version  of  the  Visibility  model. 

A.  THE  VISIBILITY  CONCEPT 

Several  models  have  been  proposed  for  predicting  the  performance  of  staring 
thermal  imaging  systems.  Among  these,  the  infrared  community  standard  model  FLIR-92 
has  gained  a  wider  acceptance.  As  mentioned  in  Chapter  III,  the  FLIR-92  model  predicts 
MRTD  below  the  Nyquist  limit  and  does  not  take  into  account  the  aliasing  effects.  It  also 
makes  a  set  of  assumptions  regarding  the  observer  eye-brain  recognition  process,  which 
brings  some  limitations  to  the  perfonnance  predictions.  Therefore  the  FLIR-92 
predictions  were  observed  to  be  overly  optimistic  at  low  spatial  frequencies  and  overly 
pessimistic  at  high  spatial  frequencies  [Ref.  17]. 

The  Visibility  model,  which  was  based  on  a  minimum  threshold  input  contrast 
parameter  and  a  contrast  reduction  factor  due  to  aliasing  and  blurring  effects,  was 
proposed  for  predicting  the  MRTD  perfonnance  of  staring  systems.  The  model  makes  no 
assumptions  about  the  observer’s  judgment  in  the  recognition  process.  Therefore  it  can 
easily  be  used  for  the  perfonnance  predictions  of  a  thermal  imaging  system  utilizing  an 
Automatic  Target  Recognizer  (ATR)  device  instead  of  a  human  observer,  which  makes 
the  Visibility  model  more  adaptable  to  the  goal  of  modeling  an  objective  MRTD 
prediction  [Ref.  3].  The  details  of  the  objective  MRTD  measurement  procedure  will  be 
discussed  later  in  Chapter  V. 
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The  Visibility  model  suggests  the  incorporation  of  a  low  frequency  limit  into  the 
MRTD  calculations  in  the  form  of 


MRTD(fz)  =MRTD(fz  ->  0)x 


1 

«(/z) 


(4.1) 


where  subscript  z  indicates  either  the  horizontal  or  the  vertical  direction,/]  is  the 
spatial  frequency  and  a (fz)  is  the  contrast  transference  parameter.  Equation  4.1  has  a 
significant  physical  appeal  because  it  is  based  on  the  fact  that  to  resolve  a  target  with  an 
input  contrast  of  zero  is  not  possible.  This  is  why  the  laboratory-based  MRTD 
measurements  extend  to  a  certain  non-zero  limit  at  low  spatial  frequencies  [Refs.  3  and 
20]. 

Depending  on  this  simple  but  important  physical  reality,  a  threshold  contrast 
value  of  AT sc  is  defined  at  the  low  frequency  limit  below  which  recognition  is  not 
possible.  ATsc  is  an  observer  and  system  dependent  term  that  inherently  comprises  all 
noise  effects  without  making  any  assumptions  about  their  natures.  For  that  reason  an 
observer  may  have  different  thresholds  for  various  thermal  imaging  systems  [Refs.  3,  19 
and  20]. 

The  block  diagram  in  Figure  4.1  shows  the  fundamental  processes  taken  into 
account  in  the  model.  First  a  four-bar  target  with  a  high  contrast  AT  is  modeled  as  a 
rectangular  function,  which  is  considered  as  the  input  to  the  system.  During  the  image 
formation  the  input  experiences  contrast  degradation  due  to  blurring  effects  at  the  optics 
while  aliasing  effects  due  to  the  sampling  process  are  introduced  in  two  dimensions  at  the 
detector  plane.  Scene  phasing  effects,  which  are  observed  due  to  the  target  misalignment 
from  the  center  of  the  detector  elements,  are  also  taken  into  account  throughout  the 
sampling  process.  Next  the  model  utilizes  the  electronics  OTF  which  brings  a  non-zero 
PTF  into  the  predictions  before  the  introduction  of  the  display  MTF  into  the  image 
reconstruction  process. 
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Figure  4.1.  The  Visibility  Model  block  diagram  [Ref.  3]. 


Finally,  the  rectangular  input  function  with  a  high  contrast  AT  appears  as  a 
smoothed-out  image  with  a  degraded  contrast  of  ATS  at  the  display  [Ref.  3].  Related  to 
the  parameters  shown  in  Figure  4.1,  the  contrast  transference  parameter  for  a  particular 
bar  spatial  frequency  is  given  as  in  Equation  4.2. 


«(/) 


AW) 

AT 


(4.2) 


Depending  on  the  physical  reasoning  of  the  Visibility  model,  the  relationship 
between  the  contrast  transference  parameter  and  MRTD  value  is  constructed  by  the 
following  equation. 
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MRTD(f  ) 


(4.3) 


^  Tsc 

«(/) 

B.  THE  CALCULATION  OL  THE  THRESHOLD  INPUT  CONTRAST 
PARAMETER  ATSC 

The  input  contrast  threshold  parameter  ATsc  in  Equation  4.3  can  be  detennined  in 
two  ways.  It  can  be  obtained  from  the  laboratory  measurements  or  it  can  be  calculated 
from  the  system  parameters  using  the  following  heuristic  fonnulation  [Ref.  19]. 


(  _2 

n  a  vh 

^  ' sc  ~  0  -j(V/v7/wO"n7/  2 

O  cr  TVH 


(Erf 


where 

SNRjhr  is  the  threshold  signal-to-noise  ratio 

Et  is  the  temporal  integration  factor  of  the  eye 

cjtvh  is  the  random  temporal  pixel  noise 

<jVh  is  the  random  spatial  noise 

Equation  4.4  apparently  shows  that  ATsc  calculations  incorporate  the  3-D  noise 
concept.  Only  temporal  pixel  noise  (aiva)  and  spatial  noise  (ctvh)  are  included  in  the 
analysis  since  they  are  the  dominant  noise  components  in  staring  systems,  which  can  be 
determined  from  measurements.  If  measurement  data  is  not  available  gtvh  can  be 
calculated  using  Equation  3.2  and  gvh  can  be  estimated  from  its  default  value  that  was 
determined  to  be  0.40gtvh  as  in  the  FLIR-92  model  [Ref.  19]. 

The  eye  temporal  integration  factor  is  also  included  as  an  enhancement  element 
into  Equation  4.4  since  the  observer  is  considered  as  a  part  of  the  system.  Therefore  this 
formula  provides  a  subjective  prediction  for  the  threshold  contrast  parameter  AT$c- 
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C.  THE  CALCULATION  OL  THE  CONTRAST  TRANSEERENCE 

PARAMETER  a(F) 

The  contrast  transference  parameter  a(f)  brings  all  frequency  dependent  system 
components  characterized  in  terms  of  the  sub-system  MTFs  into  the  MRTD  predictions. 
It  also  accounts  for  the  two-dimensional  sampling  and  aliasing  effects  that  are  inherent  in 
the  staring  thermal  imagers  [Ref.  19]. 

For  the  calculation  of  the  contrast  transference  parameter,  an  averaging  process  is 
applied  to  the  portions  of  the  output  image  corresponding  to  the  center  of  the  leftmost  bar 
and  to  the  center  of  the  leftmost  trough  lying  in  between  the  two  bars.  The  difference 
between  the  average  obtained  from  the  bar  represented  by  Tp  and  the  average  obtained 
from  the  trough  represented  by  Tt  yields  a  unaliased  output  contrast  value  as  shown  in 
Figure  4.2.  On  the  other  hand  a  set  of  alias-noise  values  is  calculated  by  subtracting  the 
aliased  output  image  from  the  unaliased  one  point  by  point.  Then  the  standard  deviation 
of  the  alias-noise  values  is  taken  to  be  the  alias-noise  parameter  denoted  by  gAlias- 
Finally  the  output  contrast  is  detennined  by  subtracting  the  alias-noise  parameter  gAlias 
from  the  unaliased  contrast  value  that  is  given  by  Equation  4.5  [Ref.  19]. 

A Ts(f  )  =  TP(f  )-TT(f  )-aAlIJS(f  ) 


As  the  input  contrast  is  set  to  one  for  all  spatial  frequencies,  the  contrast 
transference  parameter  is  calculated  as 

a{f)=ATs{f)  (4.6) 

When  this  process  is  repeated  over  the  range  of  all  spatial  frequencies,  the 
contrast  transference  parameter  is  obtained  as  a  function  of  spatial  frequency.  The  results 
obtained  in  the  previous  works  [Refs.  3  and  20]  indicated  that  the  Visibility  model 
predictions  were  reasonably  accurate.  Therefore  the  model  was  confirmed  to  be  a  simple 
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but  precise  means  of  predicting  the  performance  of  thermal  imaging  sensors  [Ref.  3]. 
This  is  the  most  significant  reason  why  the  Virtual  Thermal  Image-Processing  model  is 
based  on  the  Visibility  model. 


Figure  4.2.  A  Typical  Contrast  Calculation  Shown  on  a  Profile  Across  the 
Horizontal  Bar  Pattern  Image  Generated  at  0.65  Cycles/Mrad. 

D.  THE  LATEST  AMENDMENTS  IN  THE  VISIBILITY  MODEL 

The  Visibility  model  has  been  modified  to  enable  a  more  precise  modeling  of  the 
thermal  imaging  systems  without  changing  any  of  the  primary  concepts  of  its  origin. 
Basically,  the  modifications  can  be  discussed  under  two  main  headings:  the  modifications 
concerning  the  calculation  of  the  threshold  contrast  parameter  AT$c  and  the  modifications 
concerning  the  contrast  transference  parameter  a(f). 
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1.  Modifications  concerning  the  Calculation  of  ATsc 

For  the  calculation  of  the  threshold  contrast  parameter  ATsc,  a  correction  is 
proposed  to  Equation  4.4.  According  to  the  3-D  noise  concept  in  FLIR-92,  the  eye 
temporal  integration  factor  Ex  should  only  enhance  the  random  temporal  pixel  noise 
component  citvh  that  has  noise  in  the  temporal  direction.  Therefore  Ex  should  not  impact 
the  random  spatial  noise  (fixed  pattern  noise)  component  gvh  as  it  does  in  Equation  4.4. 
The  proposed  formulation  for  the  calculation  of  ATsc  is  given  below  in  Equation  4.7. 

2 

A  Tsc  = - SNR  THR  (cr  2tvh  Et  +  ct~vh)  (4.7) 

8 

Another  modification  affecting  ATsc  is  made  in  calculating  the  system  total 
bandwidth  A/p  used  in  the  determination  of  the  noise  component  gtvh-  The  relationship 
between  gtvh  and  A/P  was  given  as  in  Equation  3.2.  In  the  original  model,  A/P  was 
calculated  using  only  the  reconstruction  MTF.  The  amendment  suggests  including  the  eye 
MTF  into  the  A /P  calculation,  which  will  also  be  discussed  under  the  noise  analysis  topic 
of  the  Virtual  Thennal  Image-Processing  model  in  Chapter  VI. 

2.  Modifications  concerning  the  Calculation  of  a(F) 

For  the  calculation  of  the  contrast  transference  parameter  a(f),  a  new  method  has 
been  developed  in  which  a(f)  is  computed  directly  from  the  output  image.  Instead  of  an 
averaging  process,  the  maximum  values  on  each  bar  and  minimum  values  on  each  space 
are  computed  first.  Then,  contrast  values  are  calculated  by  subtracting  the  minimum  from 
the  maximum  for  each  adjacent  bar  and  space  in  the  image.  Finally,  the  minimum 
contrast  value  among  the  successive  contrast  values  in  the  image  is  named  to  be  the 
contrast  transference  value  for  that  specific  bar  spatial  frequency.  When  this  process  is 
repeated  over  the  range  of  all  spatial  frequencies,  the  contrast  transference  parameter  is 
obtained  as  a  function  of  spatial  frequency. 

Further  changes  affecting  the  contrast  transference  parameter  were  made  in  the 
calculation  of  the  sub-system  MTFs.  First,  the  methodology  of  obtaining  electronics  OTF 
was  modified.  In  the  original  model,  it  was  calculated  from  the  electronics  MTF  and 
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electronics  PTF  equations  directly.  The  new  approach  suggests  obtaining  the  electronics 
OTF  by  simply  taking  the  2-D  FFT  of  its  impulse  response  function  for  a  more  accurate 
calculation.  Second,  an  additional  detector  MTF  was  included  into  the  image 
reconstruction  process  as  treated  by  Lloyd  in  Reference  10.  The  details  for  the  image 
formation  and  image  reconstruction  are  provided  in  Chapter  6. 

All  transfer  function  equations  used  in  the  Visibility  model  predictions  will  be 
given  in  Chapter  VI  since  the  Virtual  Thermal  Image-Processing  model  utilizes  the  same 
sub-system  transfer  functions.  Further  details  of  the  revisions  accomplished  in  the 
Visibility  model  can  be  found  in  Reference  21. 

E.  THE  OBJECTIVE  MRTD  PREDICTIONS  OF  THE  VISIBILITY  MODEL 

Up  to  this  point  in  the  discussions  concerning  the  MRTD  predictions,  the  human 
observer  was  considered  as  a  part  of  the  experimental  process,  which  makes  this  type  of 
assessment  subjective  and  less  accurate  because  the  results  depend  on  the  performance  of 
the  person  making  the  measurement.  Although  experienced  observers  tend  to  produce 
repeatable  results,  it  would  be  preferred  to  make  these  measurements  such  that  the  results 
are  independent  of  the  observer  [Ref.  22]. 

With  this  motivation,  a  new  Visibility  model  is  constructed  for  predicting  the 
objective  MRTD  perfonnance  of  thermal  imaging  systems  taking  into  account  the 
enhancements  perfonned  in  the  Visibility  model  that  makes  subjective  MRTD 
predictions.  In  contrast  to  the  subjective  model,  the  objective  model  takes  the  human 
observer  and  the  display  out  of  the  loop  and  calculates  the  threshold  contrast  parameter 
AT sc  using  the  simpler  Equation  4.8. 

AT 

A^c  =  SNRmR  ^  N v  (4.8) 

In  Equation  4.8, 

SNRjhr  is  the  threshold  signal-to-noise  ratio  for  a  system  to  resolve  a  target 
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AT  is  the  temperature  difference  between  the  target  and  the 

background 

AV  is  the  differential  signal  voltage  produced  by  the  thennal  imaging 

system  due  to  the  temperature  difference  of  AT  between  the  target 
and  the  background 

Nv  is  the  noise  voltage 

Equation  4.8  is  based  on  the  assumption  that  a  differential  temperature  between 
the  target  and  the  background  would  produce  a  relative  differential  voltage  at  the  sensor, 
which  makes  the  (AT/AV)  ratio  approximately  a  constant.  The  data  obtained  from  the 
objective  MRTD  measurements  performed  in  a  laboratory  environment  also  support  this 
theory  strongly.  Another  significant  feature  of  Equation  4.8  is  that  a  SNRthr  parameter  is 
integrated  into  the  predictions  because  an  increase  in  the  SNRthr  requirement  for  a 
system  to  resolve  a  target  suggests  an  increase  in  the  system  ATsc  value,  which  makes 
Equation  4.8  applicable  for  different  system  SNRthr  requirements. 

Further  details  of  the  Visibility  model  that  makes  objective  MRTD  predictions 
can  be  found  in  Reference  21.  Both  subjective  and  objective  MRTD  predictions  of  the 
Visibility  model  will  be  compared  to  the  laboratory  measurements  in  Chapter  V. 
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V.  SUBJECTIVE  AND  OBJECTIVE  MRTD  MEASUREMENTS 
AND  THE  COMPARATIVE  ANALYSIS  OF  THE  DATA 


The  purpose  of  this  chapter  is  to  explain  the  experimental  setup  and  the 
methodology  used  in  both  subjective  and  objective  MRTD  measurements  made  in  our 
laboratory  at  the  Naval  Postgraduate  School  and  to  present  the  results.  Ultimately,  these 
results  will  also  be  compared  with  the  recently  modified  Visibility  model  and  the  FLIR- 
92  model  predictions  in  this  chapter. 

A.  THE  SUBJECTIVE  MRTD  MEASUREMENTS 

1.  Experimental  Setup 

A  staring  thermal  imaging  device,  the  Mitsubishi  Electronics  IR-M500  model, 
utilizing  a  512x512  PtSi  detector  array  was  used  in  the  subjective  MRTD  measurements. 
The  imager’s  system  parameters  are  given  in  Appendix  D. 

The  experimental  equipment  was  arranged  as  shown  in  Figure  5.1.  Composed  of  a 
front  plate  and  a  heated  back  plate,  a  target  platform  was  placed  at  a  distance  of  R  from 
the  thermal  imager.  Figure  5.2  shows  the  details  of  this  platform.  For  a  proper  alignment, 
the  imager  optical  aperture  and  the  target  front  plate  were  placed  parallel  to  each  other  at 
the  same  height.  Therefore  the  effects  of  scene  phasing  phenomenon  were  minimized  and 
a  clean  image  was  presented  to  the  observer. 

Both  plates  on  the  target  platfonn  were  made  of  shot-blasted  aluminum  to  provide 
high  thennal  conductivity  and  temperature  uniformity  on  the  surfaces  visible  to  the 
thennal  imager  and  they  were  painted  in  a  non-reflective  flat  black  to  ensure  uniform 
emissivity.  The  temperature  of  the  back  plate  was  adjusted  by  controlling  the  electrical 
current  provided  to  the  resistive  heating  elements  placed  on  its  back  surface.  On  the  other 
hand  the  front  plate  remained  near  its  ambient  temperature  since  it  was  separated  by 
10cm  from  the  back  plate.  A  thermocouple  mounted  on  each  plate  enabled  the 
temperature  difference  between  the  target  and  the  background  to  be  read  directly. 
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Figure  5.1.  The  Experimental  Setup. 


As  seen  in  Figure  5.2,  standard  four-bar  patterns  (7:1  aspect  ratio)  with  various 
sizes  were  cut  through  the  front  plates  allowing  the  observer  to  make  measurements  at 
different  spatial  frequencies.  Once  one  of  these  front  plates  was  slid  over  the  fixed 
template  in  front  of  the  back  plate,  the  thermal  imager  viewed  the  heated  back  plate 
through  the  openings  on  the  front  plate.  With  this  configuration,  the  front  plate 
represented  the  uniform  background  while  the  back  plate  represented  the  hot  target  in  the 
scene  of  interest. 
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Figure  5.2.  The  Target  Platform. 

The  measurements  were  conducted  at  an  ambient  room  temperature  of  300 
degrees  Kelvin.  The  laboratory  environment  was  organized  so  that  the  reflected  radiation 
from  other  emitting  objects  in  the  laboratory  environment  was  minimized.  Since  the 
experimental  equipment  was  static,  the  motion  effects  of  the  thermal  imager  and  the 
target  were  not  accounted  for,  which  is  consistent  with  the  static  nature  of  most  of  the 
laboratory  performance  measurements  in  the  literature. 

2.  Experimental  Methodology 

Some  preparations  were  made  before  the  actual  measurements.  The  observers 
were  trained  first  for  the  recognition  criteria.  The  target  was  called  ‘resolved’  whenever 
all  observers  fully  agreed  that  it  consisted  of  four  bars.  Therefore  the  MRTD 
measurement  process  was  based  on  a  probability  of  100%  detection,  which  corresponds 
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to  a  SNR  of  about  6.0  [Ref.  25].  Other  measurement  methodologies  utilizing  various 
detection  probabilities  suggest  different  SNR  values.  As  a  result  adjusting  the  SNR 
makes  it  easy  to  incorporate  various  measurement  methodologies  in  the  model 
predictions. 

Then  the  thermal  imager  was  focused  and  aligned  with  the  target  platform  to 
optimize  the  system  response.  After  darkening  the  laboratory,  the  observers  were 
allowed  sufficient  time  to  adapt  to  the  environmental  conditions  and  to  adjust  the  system 
gain  and  monitor  contrast  and  brightness  to  enhance  the  viewing  conditions. 

Two  sets  of  measurement  data  were  taken  during  the  experiment.  One  set  was 
taken  while  the  back  plate  was  heated  up.  The  observers  waited  until  they  could  resolve 
the  four  bars.  Then  the  second  set  of  data  was  taken  while  the  back  plate  was  left  to  cool 
down.  The  observers  waited  until  they  could  no  longer  resolve  the  four  bars.  In  both 
measurements,  the  temperature  difference  between  the  front  plate  and  the  back  plate  was 
recorded  as  the  MRTD.  These  steps  were  repeated  for  a  range  of  spatial  frequencies. 
Finally  one  set  of  MRTD  as  a  function  of  spatial  frequency  was  calculated  by  averaging 
these  data  taking  all  trial  results  into  account.  Both  the  horizontal  and  vertical  MRTD 
measurements  were  taken  repeating  the  same  process.  As  in  FLIR-92,  the  two- 
dimensional  MRTD  were  also  calculated  for  a  complete  analysis. 

3.  Measurement  Results 

The  complete  set  of  the  subjective  MRTD  measurement  results  is  given  in 
Appendix  E.  Figure  5.3  portrays  the  corresponding  MRTD  plots  as  a  function  of  spatial 
frequency.  As  mentioned  before  this  data  will  be  compared  to  the  Visibility  model  and 
the  FLIR-92  model  predictions  later  in  this  chapter. 

There  are  several  conclusions  that  should  be  pointed  out  from  the  subjective 
MRTD  measurement  results  as  follows. 

There  is  an  eye-catching  difference  between  the  horizontal  and  vertical  MRTD 
measurement  results.  System  MRTD  performance  in  the  horizontal  direction  appeared  to 


52 


be  better  than  that  in  the  vertical  direction.  This  conclusion  can  be  tied  to  the  difference 
between  the  horizontal  and  vertical  MTF  values  of  the  system. 

The  MRTDs  of  the  Mitsubishi  system  seemed  to  have  nonzero  low  spatial 
frequency  limits  that  were  measured  to  be  0.08  and  0.15  degrees  C  in  the  horizontal  and 
vertical  directions  respectively.  As  mentioned  in  Chapter  IV,  this  physical  reality  was 
used  as  the  principal  construction  milestone  of  the  Visibility  model. 

The  typical  asymptotic  characteristic  of  an  MRTD  curve  above  a  particular  spatial 
frequency  value  is  clearly  noticeable  in  the  MRTD  plots  of  the  Mitsubishi  system  as  in 
Figure  5.3.  Finally,  a  close  examination  of  the  results  and  plots  provided  in  Appendix  E 
indicates  a  significant  difference  between  the  measurements  taken  as  the  back  plate  is 
heated  up  and  the  measurements  taken  as  the  back  plate  is  cooled  down.  This  observation 
can  be  attributed  to  a  perception  experience  such  that  it  is  easier  to  track  a  bar  pattern  as  it 
disappears  after  it  has  been  resolved  than  to  recognize  a  new  pattern  as  it  starts  to  appear 
out  of  noise. 

The  tolerance  of  the  thennocouples  used  to  measure  the  plate  temperatures  was 
not  perfect.  They  had  a  measurement  precision  of  0.0556  degrees  C  that  corresponds  to 
0.1  degrees  F. 
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Figure  5.3.  Mitsubishi  System  Average  Subjective  MRTD  Plots. 

B.  THE  OBJECTIVE  MRTD  MEASUREMENTS 
1.  Experimental  Setup 

As  discussed  before  in  Chapter  IV,  the  dependence  on  human  observer  makes  the 
MRTD  measurements  subjective  and  unrepeatable.  Although  several  approaches  have 
been  proposed  in  the  literature  [Refs.  27,  28,  29  and  30],  there  is  no  standard  objective 
MRTD  measurement  procedure  that  has  been  accepted  yet.  However  the  trend  in 
operating  an  ATR  device  to  objectively  decide  the  resolvability  of  a  target  in  the  presence 
of  noise  is  becoming  an  effective  method  of  eliminating  the  human  dependence. 

The  objective  MRTD  experimental  equipment  at  the  Naval  Postgraduate  School 
was  arranged  as  shown  in  Figure  5.4.  In  contrast  to  the  subjective  MRTD  measurement 


54 


setup  the  display  and  the  observer  were  taken  out  of  the  loop  while  an  oscilloscope  was 
directly  connected  to  the  thennal  imager,  simulating  the  simplest  form  of  an  ATR  device. 
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Oscilloscope 


Uniform 

Background 


Figure  5.4.  The  Experimental  Setup. 

2.  Experimental  Methodology 

Since  ATR  based  MRTD  measurements  are  basically  related  to  SNR,  calculation 
of  noise  becomes  a  major  preliminary  process.  First  the  noise  voltage  measured  along  a 
video  scan  line  across  the  target  pattern  was  viewed  on  an  oscilloscope  screen  as  seen  in 
Figure  5.5.  Next,  the  average  noise  voltage  was  calculated  as  8.1mV  by  taking  samples 
from  the  curve  in  the  figure.  Then  the  SNR  requirement  for  the  desired  discrimination 
task  was  determined  from  the  ratio  of  the  voltages  (AVsignal/AVnoise),  where  AVnoise 
was  set  to  8.1  mV. 
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Figure  5.5.  Noise  Voltage  as  seen  on  the  Oscilloscope. 


A  typical  SNR  measurement  is  demonstrated  in  Figure  5.6.  For  example,  a  signal 
voltage  level  (AVsignal)  of  32.6mV  is  required  on  the  oscilloscope  screen  to  call  a  target 
‘resolved’  for  a  SNR  requirement  of  4.0. 

Two  sets  of  measurements  were  taken  in  the  objective  MRTD  experiment.  First, 
SNR  was  held  constant  while  the  MRTD  values  were  recorded  for  a  range  of  bar  spatial 
frequencies.  SNR  was  set  to  6.0  for  a  100%  probability  of  detection.  Second,  bar  spatial 
frequency  was  held  constant  while  the  MRTD  values  were  recorded  for  various  SNR 
requirements. 
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Figure  5.6.  A  Typical  SNR  Measurement  on  the  Oscilloscope. 

3.  Measurement  Results 

A  complete  set  of  measurement  data  is  provided  in  Appendix  E.  Figures  5.7  and 
5.8  below  show  the  plots  of  MRTD  versus  spatial  frequency  for  a  constant  SNR  of  6.0 
and  MRTD  versus  SNR  for  various  constant  spatial  frequencies  respectively. 

When  the  MRTD  values  obtained  from  the  two  sets  of  measurements  are 
compared  i.e.  for  SNR  of  6.0,  it  is  obvious  that  the  objective  MRTD  scheme  yields 
considerably  more  reproducible  results.  However,  objective  MRTD  values  turned  out  to 
be  much  higher  than  the  MRTD  values  obtained  from  the  subjective  measurements  i.e. 
for  SNR  value  of  6.0.  This  conclusion  suggests  that  the  eye  is  an  outstanding  detector  of 
signal  out  of  noise.  The  objective  measurement  results  will  be  compared  to  the  objective 
predictions  of  the  Visibility  model  also  in  this  chapter. 
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Figure  5.7.  Objective  MRTD  vs.  Spatial  Frequency  for  a  SNR  of  6.0. 
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Figure  5.8.  Objective  MRTD  vs.  SNR  for  Constant  Spatial  Frequencies. 


C.  COMPARATIVE  DATA  ANALYSIS 

The  purpose  of  this  subsection  is  to  make  a  comparative  analysis  of  the  Visibility 
model  predictions  to  the  measured  data.  Furthermore  this  analysis  will  validate  the 
applicability  of  the  Visibility  concept  to  staring  thermal  imagers. 

First,  the  subjective  MRTD  measurement  results  will  be  compared  to  both  the 
Visibility  model  and  the  FLIR-92  model  subjective  MRTD  predictions.  Second,  the 
objective  MRTD  measurement  results  will  be  compared  to  the  Visibility  model  objective 
predictions.  All  prediction  results  of  both  models  are  presented  graphically  in  Appendix 
F. 
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1.  Comparison  of  the  Subjective  MRTD  Measurements  to  the 
Predictions 

In  Figures  5.9,  5.10  and  5.11,  horizontal,  vertical  and  the  2-D  MRTD  predictions 
of  the  Visibility  and  the  FLIR-92  models  are  compared  to  the  MRTD  measurement 
results.  The  2-D  MRTD  measurement  data  is  calculated  by  geometrical  averaging  of  the 
measured  horizontal  and  vertical  MRTD  values,  in  the  same  manner  that  the  other  2-D 
MRTD  values  were  obtained. 

As  seen  in  Figure  5.9,  the  horizontal  MRTD  predictions  of  both  models  are 
observed  to  be  optimistic  at  low  spatial  frequencies  and  pessimistic  at  high  spatial 
frequencies.  On  the  other  hand,  both  models  seem  to  give  optimistic  results  for  the 
vertical  and  the  2-D  MRTD  values  as  demonstrated  in  Figures  5.10  and  5.11  respectively. 


Figure  5.9.  Horizontal  MRTD  Comparison. 
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Figure  5.10.  Vertical  MRTD  Comparison. 

In  Figures  5.9  and  5. 10,  it  is  noticeable  that  the  Visibility  and  the  FLIR-92  models 
tend  to  give  quite  close  predictions  in  the  horizontal  and  vertical  directions  up  to  the 
Nyquist  frequency  limit.  It  is  also  observable  that  the  FLIR-92  predictions  are  limited  to 
the  Nyquist  frequency  in  the  corresponding  direction.  However  the  Visibility  model 
predictions  extend  beyond  this  point,  because  the  actual  laboratory  measurements  suggest 
that  resolution  is  still  possible. 
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Figure  5.11.  2-D  MRTD  Comparison. 

Figures  5.12,  5.13  and  5.14  are  given  to  show  the  tendency  in  MRTD  predictions 
in  the  low  spatial  frequency  limit.  It  is  clearly  observable  that  the  Visibility  model 
calculates  the  low  spatial  frequency  limit  that  was  originally  defined  as  ATsc  (threshold 
input  contrast)  in  the  model  quite  precisely.  However  the  FLIR-92  predictions  tend  to  go 
to  zero  in  the  low  spatial  frequency  region. 
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Figure  5.12.  Horizontal  MRTD  Comparison  at  Low  Spatial  Frequencies. 
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Figure  5.13.  Vertical  MRTD  Comparison  at  Low  Spatial  Frequencies. 
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Figure  5.14.  2-D  MRTD  Comparison  at  Low  Spatial  Frequencies. 

Finally,  it  can  be  said  that  the  Visibility  model  predictions  are  in  better  agreement 
with  the  measured  results  than  are  the  FLIR-92  model  predictions,  as  noticed  in  the 
figures  given  in  this  section. 

2.  Comparison  of  the  Objective  MRTD  Measurements  to  the  Predictions 

Since  there  is  no  objective  MRTD  predictor  version  of  the  FLIR-92  model  yet,  the 
measurement  results  are  only  compared  to  the  Visibility  model  predictions.  First,  the 
measurement  results  are  compared  to  the  predictions  when  the  SNR  is  held  constant 
(SNR=6.0)  as  also  demonstrated  in  Figure  5.15.  Although  the  Visibility  model 
predictions  seem  to  be  more  optimistic  than  the  measured  results,  there  is  a  reasonable 
agreement  between  the  two  curves.  Moreover  the  measurement  and  the  prediction 
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nonzero  low  frequency  limits  are  quite  close  to  each  other,  indicating  that  the  basic 
physical  reasoning  of  the  Visibility  model  is  correctly  constructed. 


Figure  5.15.  Objective  MRTD  Comparison  (Constant  SNR  =  6.0). 

Second,  the  measurement  results  are  compared  to  the  predictions  as  a  function  of 
SNR  when  bar  spatial  frequency  is  held  constant  as  demonstrated  in  Figures  5.16,  5.17, 
5.18  and  5.19.  The  spatial  frequencies  were  chosen  to  be  0.12,  0.24,  0.46  and  0.56 
cycles/mrad.  Each  figure  clearly  portrays  the  agreement  between  the  measured  and  the 
predicted  results.  However  the  best  agreement  is  observed  at  0.56  cycles/mrad  bar  spatial 
frequency  (Figure  5.19),  which  was  also  the  maximum  bar  spatial  frequency  utilized  in 
these  measurements/predictions  (Figure  5.15). 


66 


Figure  5.16.  Objective  MRTD  vs.  SNR  Comparison  for  Constant  Spatial 

Frequency  (0.12  cycles/mrad). 
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Figure  5.17.  Objective  MRTD  vs.  SNR  Comparison  for  Constant  Spatial 

Frequency  (0.24  cycles/mrad). 
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Figure  5.18.  Objective  MRTD  vs.  SNR  Comparison  for  Constant  Spatial 

Frequency  (0.46  cycles/mrad). 
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Figure  5.19.  Objective  MRTD  vs.  SNR  Comparison  for  Constant  Spatial 

Frequency  (0.56  cycles/mrad). 


This  analysis,  although  limited  in  scope,  shows  that  the  Visibility  model  can 
provide  reasonably  accurate  performance  predictions  of  thermal  imaging  systems. 
Therefore  the  laboratory  confirmation  of  the  model  for  both  subjective  and  objective 
MRTD  schemes  is  considered  to  be  successfully  accomplished. 
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VI.  VIRTUAL  THERMAL  IMAGE-PROCESSING  MODEL 


The  purpose  of  this  chapter  is  to  introduce  the  basics  of  a  new  virtual  thermal 
image-processing  model,  which  is  proposed  in  this  thesis.  This  visualization  tool  is  based 
on  an  earlier  simulation  program,  the  Visibility  model. 

A.  INTRODUCTION 

Recently,  the  U.S.  Army  Night  Vision  Electronics  Sensors  Directorate  has  begun 
performing  Virtual  MRTD  experiments  pursuing  a  concept  known  as  simulation  based 
modeling.  As  a  result,  a  simulation  program  based  on  the  FLIR-92  model  was  developed. 
In  the  experiments,  a  scanning  thermal  imaging  system  operating  in  the  mid-wave 
infrared  band  was  simulated  and  the  effects  of  blurring  and  noise  on  the  image  were 
reported  [Ref.  23], 

Inspired  by  this  preceding  work,  a  new  Virtual  Image-Processing  Model  that  is 
based  on  the  Visibility  model  has  been  developed  at  the  Naval  Postgraduate  School.  Like 
the  Visibility  model,  the  Virtual  Thermal  Image-Processing  model  is  written  using  the 
MATLAB  computational  software.  It  takes  the  current  thermal  imaging  system-modeling 
concerns  such  as  noise,  scene  phasing,  sampling  and  aliasing  into  account  and  creates 
visual  images  comparable  to  those  that  can  be  obtained  with  the  actual  system  being 
modeled.  This  allows  the  user  to  analyze  and  evaluate  separately  the  effects  of  noise, 
scene  phasing,  sampling  and  aliasing  on  the  final  imagery.  Moreover  the  model  provides 
control  over  the  system  parameters  for  the  analyst,  which  might  be  helpful  in  perfonning 
realistic  tests  in  the  design  of  the  systems  under  evaluation. 

B.  THE  LOGICAL  FLOW  OF  THE  MODEL 

At  the  front  end  of  the  Virtual  Thennal  Image-Processing  model,  the  system 
parameters  of  a  particular  thermal  imaging  device  must  initially  be  loaded.  For  this  thesis, 
the  parameters  of  the  Mitsubishi  Electronics  IR-M500  Thermal  Imager,  which  is  also 
used  in  the  laboratory  measurements,  are  entered  as  previously  given  in  Appendix  D. 

To  provide  a  more  convenient  implementation  of  the  effects  of  relevant  concerns 
such  as  sampling  and  aliasing  in  the  thermal  imaging  systems,  the  program  primarily 
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operates  in  the  two-dimensional  spatial  frequency  domain.  Hence  the  two  dimensional 
Fast  Fourier  Transform  (2-D  FFT)  of  the  bar  pattern  is  taken  to  obtain  a  four-bar  target 
spectrum. 
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Figure  6.1.  The  Input  Bar  Pattern  is  Represented  by  a  256  x  256  Matrix  of  Values 

Zero  and  One. 


The  Virtual  Thermal  Image-Processing  Model  starts  with  generating  a  numerical 
model  of  the  input  four-bar  pattern  by  256  x  256  points  of  values  zero  and  one  as  shown 
in  Figure  6.1.  As  the  images  at  various  spatial  frequencies  are  analyzed,  the  number  of 
points  (256x256)  representing  the  input  bar  pattern  stays  constant,  but  the  number  of 
points  per  unit  angular  subtense  (mrad)  is  varied.  Therefore  the  axes  on  the  output  images 
are  properly  scaled  with  respect  to  the  desired  bar  spatial  frequency.  Then  the  model 
operates  for  one  spatial  frequency  at  a  time  and  loops  through  all  spatial  frequencies  of 
interest  to  provide  an  output  image  for  each  spatial  frequency.  It  works  for  both 
horizontal  and  vertical  bar  patterns  simultaneously. 
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Figure  6.2.  Block  Diagram  of  the  Virtual  Thermal  Image-Processing  Model. 
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In  the  next  stage,  the  model  computes  an  Image  Fonnation  MTF  that  is  multiplied 
with  the  four-bar  target  spectrum  in  the  spatial  frequency  domain  as  shown  in  the  block 
diagram  in  Figure  6.2.  Then  white  noise  is  generated  to  account  for  the  overall  system 
noise.  After  taking  its  2-D  FFT,  the  noise  spectrum  is  obtained  and  simply  added  to  the 
filtered  target  spectrum.  Further  discussions  about  the  calculation  and  the  generation  of 
white  noise  will  be  given  in  the  noise  analysis  of  this  chapter. 

The  sampling  process  is  applied  to  the  pre-filtered,  noise  added  target  spectrum  in 
the  spatial  frequency  domain  as  performed  in  the  Visibility  model.  The  net  result  of  this 
process  is  the  occurrence  of  the  sampled  spectrum  replicas  at  the  multiples  of  the 
sampling  frequency  in  both  horizontal  and  vertical  directions. 

Next  the  Virtual  Thermal  Image-Processing  Model  calculates  the  Image 
Reconstruction  filter  frequency  response  and  then  multiplies  it  with  the  aliased  target 
spectrum  in  the  spatial  frequency  domain.  Finally  taking  the  inverse  2-D  FFT  of  the 
system  output  spectrum  gives  a  visual  four-bar  target  representation  at  that  particular 
spatial  frequency. 

C.  ANALYSIS  OF  THE  SUB-SYSTEM  MTFS 

Since  the  Virtual  Thermal  Image-Processing  model  is  based  on  the  Visibility 
model,  it  utilizes  the  same  sub-system  MTFs  as  mentioned  before.  The  model  calculates 
the  overall  system  response  by  simply  multiplying  them  in  the  frequency  domain. 


As  shown  in  Figure  6.2,  the  sub-system  MTFs  are  categorized  into  two  main 
groups,  the  Image  Formation  MTF  and  Image  Reconstruction  MTF.  The  Image 
Formation  MTF  includes  the  diffraction-limited  optical  MTF,  the  geometric  blurred 
optical  MTF,  the  detector  spatial  MTF  and  an  optional  scene  phasing  MTF,  which  are 
given  in  the  following  equations. 
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(6.3) 
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Where 

f  is  the  two-dimensional  spatial  frequency  in  cy/mrad 

fopT  is  the  optical  cutoff  spatial  frequency  in  cy/mrad 

ctblur  is  the  standard  deviation  of  the  blur  spot  diameter  (ctblur=0  in  this 
study) 

z  is  the  direction  of  interest  either  horizontal  or  vertical 

8  is  the  detector  angular  subtense  in  mrad 

9  is  the  scene  phase  angle  in  rad  (0=0  in  this  study) 

fs  is  the  detector  sampling  frequency  in  cy/mrad 

Since  the  system  optics  is  assumed  to  have  a  circular  geometry,  the  diffraction 
limited  and  the  geometric  blurred  optical  MTFs  are  symmetric  in  two  dimensions.  The 
detector  MTF  might  also  be  symmetric  depending  on  its  geometry. 

Finally,  the  net  Image  Formation  transfer  function  is  calculated  as: 
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The  plot  of  the  Image  Formation  MTF  for  the  imaging  system  being  modeled  is 
given  below  in  Figure  6.3. 
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Figure  6.3.  Image  Formation  MTF  of  the  Modeled  Imaging  System. 

The  second  group  of  sub-system  MTFs  is  called  the  Image  Reconstruction  MTF, 
which  is  obtained  from  the  cascaded  form  of  the  electronics  OTF,  the  detector  MTF  and 
the  display  MTF.  The  electronics  OTF  is  calculated  by  taking  the  2-D  FFT  of  its  impulse 
response  function,  which  is  given  in  Equation  6.6,  as  discussed  previously  under  the  topic 
of  the  latest  amendments  in  the  Visibility  model.  The  additional  detector  MTF  has  the 
same  formulation  as  given  in  the  image  formation  process,  as  suggested  by  Lloyd.  The 
equations  used  in  the  calculation  of  the  Image  Reconstruction  MTF  are  given  as: 


76 


(6.6) 
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Where 

fc  is  the  electronics  cutoff  frequency  in  cy/mrad 

Z  is  the  angular  space  in  mrad 

8  is  the  detector  angular  subtense  in  mrad 

f  is  the  two  dimensional  spatial  frequency  in  cy/mrad 

a  is  the  standard  deviation  of  the  monitor  Gaussian  blur  spot  diameter 

It  is  important  to  note  that  with  the  calculation  of  the  electronics  OTF  from  its 
impulse  response,  the  non-zero  PTF  of  the  electronics  is  also  included  into  the 
calculations.  The  Image  Reconstruction  transfer  function  can  be  written  in  the  cascaded 
form  of: 
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The  plot  of  the  Image  Reconstruction  MTF  for  the  imaging  system  being  modeled 
is  given  in  Figure  6.4. 
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Figure  6.4.  Image  Reconstruction  MTF  of  the  Modeled  Imaging  System. 

Finally,  the  overall  system  MTF  can  be  calculated  from  Equation  6.10.  Figure  6.5 
given  below  shows  the  plot  of  the  overall  system  MTF. 
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Figure  6.5.  Overall  System  MTF  of  the  Modeled  Imaging  System. 

D.  NOISE  ANALYSIS 

The  Virtual  Thermal  Image-Processing  model  incorporates  the  updated  noise 
treatment,  the  3-D  noise  concept,  in  the  calculations  as  the  Visibility  model  does.  In  both 
models,  only  the  temporal  pixel  noise  component  Otvh  and  the  spatial  noise  component 
Gvh  are  included  in  the  analysis  since  they  are  the  dominant  noise  elements  in  the  staring 
thermal  imaging  systems.  As  discussed  previously  in  Chapter  III  gtvh  is  the  only  noise 
component  that  can  be  determined  either  from  measurements  or  from  system  parameters, 
while  the  rest  of  the  noise  components  can  only  be  detennined  from  measurements.  If 
measurement  data  is  not  available  Otvh  can  be  calculated  using  Equation  3.2  on  page  46. 
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The  most  significant  feature  of  Equation  3.2  is  that  it  utilizes  the  actual  total 
bandwidth  A fP  of  the  system  being  modeled.  The  Virtual  Thermal  Image-Processing 
model  calculates  this  quantity  using  the  following  equation  as  the  Visibility  model  does. 


ou 

A/p  =  J  |  H 


IMA  GERECONSTR  UCTION 


(fWErE(f)  df 


(6.11) 


Equation  6.11  is  the  only  place  where  a  calculation  concerning  the  observer  is 
included  into  the  model  predictions.  After  the  detennination  of  the  aTVH  parameter  from 
Equation  3.2,  the  overall  system  noise  can  be  calculated  using  Equation  3.1.  Since  the 
noise  components  other  than  Gtvh  and  <jvh  are  neglected  and  Gvh  is  estimated  from  its 
default  value  that  is  set  at  OAOotvh  as  in  the  FLIR-92  model,  the  overall  system  noise  can 
be  calculated  from  Equation  6. 12  given  below. 

® TOTAL  ~  O TVH  (l  +  0. 1  6)  (6. 12) 


Once  ototal  is  calculated,  the  model  generates  Gaussian  white  noise  with  zero 

2 

mean  and  a  variance  of  ct  white  based  on  the  relationship  given  in  the  following  equation 
below. 


a 


WHITE 


(  f  ^ 

J  MAX 

\^FX^FY  J 


a 


TOTAL 


(6.13) 


Where 


/ max  is  the  maximum  spatial  frequency  that  can  be  represented  using  the 

2-D  FFT  in  cy/mrad 

Bpx  is  the  noise  equivalent  bandwidth  in  the  x  direction 

Bfy  is  the  noise  equivalent  bandwidth  in  the  y  direction 
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As  given  in  Table  3.2  previously,  the  detector  is  the  main  source  of  the  noise 
components  cttvh  and  avH  in  staring  thermal  imaging  systems.  Therefore,  the  white  noise 
spectrum  is  inserted  into  the  model  at  the  detector  to  visualize  the  effect  of  overall  system 
noise  in  the  output  images. 

On  the  other  hand,  the  MATLAB  software  brings  a  contrast  limitation  such  that 
the  image  contrast  values  should  be  in  the  range  of  zero  to  one  for  a  proper  visualization. 
Otherwise  the  output  images  of  the  model  may  mislead  the  analyst  to  erroneous 
conclusions.  Therefore  the  amplitude  of  the  points  forming  an  image  should  be  kept  in 
the  range  of  zero  to  one.  However  addition  of  Gaussian  white  noise  at  the  detector 
produces  images  such  that  the  condition  stated  above  is  violated. 

A  practical  solution  to  this  problem  suggests  using  the  function  given  in  Equation 
6.14,  which  appropriately  nonnalizes  the  amplitude  of  all  points  in  the  image  such  that 
they  are  in  the  range  of  interest. 


F(x)  =  0.5 


tan  1  (^(jc  —  0.5))  + 1 


(6.14) 


Where 

x  is  the  input  variable 

F(x)  is  the  function  output 

a  is  a  scaling  constant  set  at  5. 

To  provide  a  better  insight,  the  plot  of  this  function  is  given  in  Figure  6.6.  The 
horizontal  axis  shows  the  input  values  while  the  vertical  axis  shows  the  corresponding 
output  values.  Although  the  input  variable  x  seems  to  be  limited  in  the  range  from  —1.5 
up  to  +1.5,  it  takes  values  from  -00  up  to  +00  yielding  output  values  between  0  and  1 
for  a  a  value  of  5. 

The  effect  of  noise  in  the  output  imagery  gets  diminished  due  to  the  application  of 
the  normalization  function  later  in  the  model.  To  compensate  for  this  side  effect,  the 


81 


amplitude  of  the  noise  spectrum  has  to  be  amplified  prior  to  its  addition.  Therefore  the 
derivative  of  F(x)  is  calculated  either  at  x=0  or  at  x=\  as  function  of  a,  as  given  in 
Equation  6.15. 


dF(x ) 


dx 


x— »0orl 


a 

1  +  0.25  a2 


(6.15) 


The  transformation  will  reduce  the  noise  by  a  factor  predicted  from  Equation 
6.15.  Finally,  the  l/g(a)  value  is  multiplied  with  the  noise  spectrum  as  an  amplification 
factor  which  was  calculated  to  be  1.45  for  a=5. 


Figure  6.6.  Normalization  Function  that  Keeps  the  Amplitude  of  Points  in  the 

Image  Between  Zero  and  One. 
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E.  SAMPLE  IMAGES 

The  purpose  of  this  part  is  to  illustrate  some  sample  images  using  the  Virtual 
Thermal  Image-Processing  model.  The  MATLAB  code  of  the  model  used  to  create  all 
images  presented  in  this  thesis  is  given  in  Appendix  G. 

Figure  6.7  and  6.8  show  the  horizontal  and  the  vertical  four-bar  target  patterns 
created  as  inputs  to  the  modeled  thermal  imager. 

In  Figures  6.9  and  6.10,  the  output  images  generated  for  0.45  cy/mrad  bar  spatial 
frequency  from  the  input  patterns  given  in  Figures  6.7  and  6.8  are  presented  respectively. 
Since  the  input  patterns  are  generated  only  once  independent  of  spatial  frequency,  they 
remain  the  same.  But  all  sub-system  MTFs,  the  sampling  process  and  the  image  axes 
throughout  the  model  are  scaled  appropriately  for  each  spatial  frequency.  Therefore  the 
output  images  are  created  with  respect  to  the  target  bar  pattern  sizes  as  seen  in  the 
following  figures. 

A  further  analysis  of  the  thennal  images  obtained  from  the  Virtual  Thennal 
Image-Processing  model  will  be  made  in  Chapter  VII. 
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Figure  6.7. 


The  Horizontal  Four-Bar  Pattern  Created  as  an  Input  to  the  Model. 


Figure  6.8.  The  Vertical  Four-Bar  Pattern  Created  as  an  Input  to  the  Model. 
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Figure  6.9.  The  Horizontal  Output  Image  Created  for  0.45  cy/mrad. 


x  space,  mrad 


Figure  6.10.  The  Vertical  Output  Image  Created  for  0.45  cy/mrad. 
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VII.  THE  ANALYSIS  OF  THE  VIRTUAL  THERMAL  IMAGE- 
PROCESSING  MODEL  OUTPUTS  AND  CONCLUSIONS 

The  purpose  of  this  chapter  is  to  analyze  the  output  images  of  the  four-bar  targets 
produced  by  the  Virtual  Thermal  Image-Processing  model  and  to  provide  useful  insights 
into  some  of  the  problems  encountered  in  modeling  the  TISs  such  as  sampling,  aliasing 
and  noise.  Moreover  the  effects  of  electronics  PTF  will  be  examined  in  the  imagery. 

The  analysis  will  be  made  under  these  sub-topics: 

•  Aliasing  effects  below  the  Nyquist  limit 

•  Aliasing  effects  above  the  Nyquist  limit 

•  Exploration  of  distortion  due  to  aliasing 

•  Electronics  PTF  effects 

Finally,  this  chapter  will  close  with  some  conclusions  and  suggested  directions  for 
future  studies. 

A.  ALIASING  EFFECTS  BELOW  THE  NYQUIST  LIMIT 

As  discussed  in  Chapter  II,  sampling  theory  suggests  that  a  frequency  component 
in  the  original  image  above  the  Nyquist  limit  appears  as  a  lower  frequency  component 
after  sampling  as  a  result  of  the  phenomenon  called  ‘aliasing’.  However,  experiments 
with  the  Virtual  Thermal  Image-Processing  model  showed  that  aliasing  effects  were  also 
present  in  the  image  of  four-bar  patterns  at  spatial  frequencies  below  the  Nyquist  limit. 
Figures  7.1  and  7.2  demonstrate  the  horizontal  four-bar  target  images  at  0.65  cycles/mrad 
including  and  excluding  the  effects  of  aliasing  respectively. 

The  spatial  sampling  frequencies  in  the  horizontal  and  vertical  directions  for  the 
thennal  imaging  system  being  modeled  were  calculated  as  1.92  and  2.5  cycles/mrad 
respectively.  Thus,  aliasing  effects  were  not  expected  to  appear  below  0.96  cycles/mrad 
bar  spatial  frequency.  However  the  images  in  Figures  7.1  and  7.2  revealed  the  aliasing 
effects  below  the  Nyquist  limit. 
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Figure  7.1.  Modeled  Four-Bar  Target  Image  with  Aliasing  At  0.65  Cycles/Mrad. 
Note  the  Apparent  Broadening  of  the  Slit  Images  and  Appearance  of  Local  Minima. 


x  space,  mrad 


Figure  7.2.  Four-Bar  Target  Image  without  Aliasing  at  0.65  Cycles/Mrad.  Note 
that  there  is  No  Such  Broadening  Effect  in  this  Image. 
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An  analytic  methodology  employing  a  spectrum  analysis  was  developed  to 
support  the  observed  effects  of  aliasing  below  the  Nyquist  limit  as  presented  above  [Ref. 
31].  Therefore  Equation  7.1,  which  gives  the  two-dimensional  Fourier  transform  of  the 
four-bar  pattern,  was  used  in  generating  the  horizontal  bar  spectrum.  Reference  31 
illustrates  Equation  7.1  in  a  different  form  after  letting  kx=2nfx  and  ky=2nfy. 


r  1  A 


V*fxJ 
^  1  ^ 
J 


(-  sin  (nfxW  )+  sin  (3 nfxW  )-  sin  (5 nfxW  )+  sin  (infJV  )) 


sin 


(l7TfrW  ) 


(7.1) 


Where 

/x  is  the  spatial  frequency  along  the  horizontal  direction  in  cycles/mrad 

fy  is  the  spatial  frequency  along  the  vertical  direction  in  cycles/mrad 

W  is  the  angular  width  of  a  bar  in  mrad 

A  plot  of  the  horizontal  four-bar  spectrum  along  the  x-axis  (/y=0)  at  0.65 
cycles/mrad  is  given  in  Figure  7.3.  The  original  spectrum,  which  is  centered  at  zero 
cycles/mrad,  repeats  itself  at  intervals  of  the  sampling  frequency.  More  importantly,  the 
figure  reveals  the  fact  that  the  replicas  overlap  with  the  original  target  spectrum  although 
the  four-bar  target  frequency  is  below  the  Nyquist  limit.  The  result  of  this  analysis  is  in 
agreement  with  the  previous  result  obtained  from  the  Virtual  Thennal  Image-Processing 
model.  Therefore  it  is  important  to  note  that  aliasing  effects  should  be  considered  for  bar- 
targets  at  four-bar  spatial  frequencies  also  below  the  Nyquist  limit. 
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Figure  7.3.  Horizontal  Four-Bar  Spectrum  Along  X-Axis  (/y=0)  at  0.65 
Cycles/Mrad  Showing  the  Overlap  of  the  Aliased  Spectra. 


A  last  observation  in  the  region  below  the  Nyquist  frequency  concerns  the 
contrast  level  in  the  target  images  generated  by  the  model.  Figure  7.4  demonstrates  the 
plot  of  a  profile  across  the  four-bar  target  and  image  centers  for  the  0.65 -eye les/mrad 
target.  Although  there  is  a  general  interpretation  that  aliasing  causes  a  visual  degradation 
effect  in  the  image,  this  analysis  shows  that  aliasing  may  actually  cause  noticeable 
contrast  enhancements  in  the  four-bar  target  imagery  for  the  bar-target  spatial  frequencies 
below  the  Nyquist  limit. 
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Figure  7.4.  A  Profile  across  the  Target  and  the  Image  Centers  at  0.65 
Cycles/Mrad.  Note  the  Enhancement  of  Contrast  (between  8%-20%)  due  to 

Aliasing. 

B.  ALIASING  EFFECTS  ABOVE  THE  NYQUIST  LIMIT 

To  account  for  aliasing  in 

modeling  thermal  imagers,  a  number  of  approaches  have  been  proposed  up  to 
now.  S.K.  Park  and  R.  Hazra  accomplished  one  of  the  notable  works  in  that  area 
presenting  a  quantitative  and  a  qualitative  assessment  of  aliasing  as  noise.  They  made  a 
theoretical  and  an  experimental  argument  suggesting  that  aliasing  should  be  treated  as 
signal-dependant,  additive  noise  [Ref.  26]. 
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Figure  7.5.  Four-Bar  Target  Image  with  Noise  and  Aliasing  At  0.05  Cycles/Mrad: 

Noise  is  Dominant  at  Low  Spatial  Frequency. 

However,  observations  made  in  the  experiments  with  the  NPS  Virtual  Thermal 
Image-Processing  model  show  that  aliasing  effects  manifest  themselves  in  much  different 
manner  from  noise  in  the  generated  images.  Figure  7.5  above  demonstrates  the  effects  of 
noise  that  were  observed  to  be  dominant  at  low  four-bar  spatial  frequencies  (0.05 
cycles/mrad).  The  tendency  of  masking  the  four-bar  pattern  was  the  most  notable 
characteristic  of  noise  as  also  observed  from  the  figure. 

On  the  other  hand,  aliasing  effects,  which  become  more  severe  at  four-bar  spatial 
frequencies  above  the  Nyquist  limit,  created  distortions  in  the  imagery.  Figure  7.6  shows 
a  significant  distortion  in  the  imagery  for  a  four-bar  spatial  frequency  of  1.10 
cycles/mrad.  Therefore  the  four-bar  pattern  appears  as  a  three-bar  pattern.  However, 
noise  alone  does  not  distort  the  image  at  the  same  four-bar  spatial  frequency;  all  bars  are 
resolvable  as  given  below  in  Figure  7.7. 
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Figure  7.6.  Four-Bar  Target  Image  with  Noise  and  Aliasing  At  1.10  Cycles/Mrad: 
Aliasing  is  Dominant  above  the  Nyquist  Frequency:  Note  the  Distortion. 


x  space,  mrad 


Figure  7.7.  Four-Bar  Target  Image  with  Only  Noise  (No  Aliasing)  at  1.10 
Cycles/Mrad:  Noise  Alone  Does  Not  Distort  the  Image. 
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This  analysis  indicates  that  aliasing  effects  cannot  be  adequately  represented  as 
noise  alone. 

C.  EXPLORATION  OF  DISTORTION  DUE  TO  ALIASING 

Experiments  with  the  Virtual  Thermal  Image-Processing  model  showed  that  there 
is  a  noticeable  difference  in  terms  of  distortion  due  to  aliasing  between  the  horizontal  and 
the  vertical  bar  pattern  images  for  the  same  spatial  frequency.  Figures  7.8  and  7.9  depict 
the  images  of  the  horizontal  and  the  vertical  bar  patterns  at  1.15  cycles/mrad  respectively 
(Nyquist  frequency  is  0.96  cycles/mrad).  The  vertical  bar  pattern  image  seems  to  be 
much  more  distorted  due  to  aliasing  when  compared  to  the  horizontal  bar  pattern  image. 

Further  analysis  showed  that  different  sampling  frequencies  in  the  horizontal  and 
vertical  directions  lead  to  the  differences  in  the  level  of  distortion  in  the  corresponding 
directions.  From  the  Mitsubishi  system  parameters,  the  sampling  frequencies  were 
calculated  as  1.92  and  2.5  cycles/mrad  in  the  horizontal  and  the  vertical  directions 
respectively.  Therefore  the  replicas  of  the  target  spectrum  are  expected  to  overlap  with 
the  original  target  spectrum  to  a  greater  extend  in  the  horizontal  direction  because  the 
sampling  frequency  in  that  direction  is  much  smaller.  At  this  point,  it  will  be  more 
appropriate  to  make  a  spectrum  analysis  to  explain  this  phenomenon. 
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Figure  7.8.  Horizontal  Four-Bar  Pattern  at  1.15  Cycles/Mrad  (above  Nyquist): 

Including  Aliasing. 


x  space,  mrad 


Figure  7.9.  Vertical  Four-Bar  Pattern  at  1.15  Cycles/Mrad:  Clear  and  Even 

Pattern  of  Horizontal  Bars. 
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The  horizontal  four-bar  spectra  along  the  x  and  y  directions  are  given  in  Figures 
7.10  and  7.11  respectively.  In  both  figures,  the  original  target  spectrum  centered  at  the 
origin  between  the  replicas  on  two  sides,  demonstrate  the  degree  of  overlap  due  to 
aliasing.  It  is  significant  that  the  overlap  of  the  spectra  along  the  x-axis  is  much  more 
when  compared  to  the  overlap  of  the  spectra  along  the  y-axis.  Since  the  sampling 
frequency  in  each  direction  is  independent  of  target  pattern  rotation,  stronger  aliasing 
effects  still  apply  along  the  x-axis  for  the  vertical  bar  pattern. 
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Figure  7.10.  Horizontal  Four-Bar  Spectrum  Along  X-Axis  at  1.15  Cycles/Mrad. 
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Hence,  strong  aliasing  effects  in  the  horizontal  direction  make  the  two  bars  in  the 
middle  appear  to  merge  as  in  Figure  7.8,  while  causing  the  vertical  bar  pattern  appear  to 
be  perpendicularly  dissected  as  in  Figure  7.9. 
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Figure  7.11.  Horizontal  Four-Bar  Spectrum  Along  Y-Axis  at  1.15  Cycles/Mrad. 


The  experiments  with  the  Virtual  Thermal  Image-Processing  model  also  showed 
the  presence  of  an  aliasing  pattern  that  manifests  itself  as  ripples  in  the  produced  images. 
Figure  7.12  below  obviously  demonstrates  this  horizontal  ripple  effect  brought  about  by 
aliasing  in  a  vertical  four-bar  image. 
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Figure  7.12.  Vertical  Four-Bar  Image  at  0.70  Cycles/Mrad. 

In  a  further  analysis  the  number  of  ripples  present  in  the  patterns  brought  about  by 
aliasing  was  counted  for  all  spatial  frequencies  in  the  vertical  four-bar  images.  Since  the 
image-reconstruction  detector  MTF  diminishes  the  ripple  effects  from  the  horizontal  bar 
pattern  images,  only  the  vertical  bar  pattern  images  are  used  for  analysis  purposes.  A 
detailed  discussion  of  the  image  reconstruction  detector  MTF  filtering  effects  is  provided 
in  Appendix  H. 

The  data  obtained  from  the  ripple  counting  process  in  the  aliasing  pattern  is 
presented  in  a  tabular  form  in  Appendix  I.  Figure  7.13  demonstrates  two  different  curves 
plotted  from  this  data.  The  upper  curve  in  the  figure  illustrates  the  number  of  ripples 
versus  bar  spatial  frequency  relationship,  while  the  lower  curve  confirms  that  the  aliasing 
pattern  stays  almost  steady  independent  of  the  bar  spatial  frequency  in  all  produced 
images. 
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Figure  7.13.  Aliasing  Pattern  Analysis  Results. 


At  this  point  a  similar  analysis  showed  that  there  is  a  reasonable  correlation 
between  the  modeled  system’s  (Mitsubishi)  sampling  frequency  and  the  number  of 
ripples  in  the  aliasing  pattern  images.  Figure  7.14  graphically  illustrates  the  data  obtained 
from  this  analysis.  The  complete  set  of  data  is  provided  in  Appendix  J.  This  figure 
verifies  how  close  the  sampling  frequency  can  be  derived  from  the  aliasing  pattern 
present  in  all  virtual  images  utilizing  a  ripple  analysis.  In  other  words,  the  ripple  effect 
brought  by  aliasing  into  the  thermal  images  is  a  function  of  the  spatial  sampling 
frequency  due  to  the  detector  array. 
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Figure  7.14.  Comparison  of  the  Sampling  Frequency  Results  Obtained  from  the 
Analysis  to  the  Sampling  Frequency  of  the  Modeled  System. 
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D.  ELECTRONICS  PTF  EFFECTS 

In  this  subsection  the  electronics  PTF  effects  will  be  examined  and  compared  to 
the  aliasing  effects  in  the  generated  images.  Figures  7.15  and  7.16  show  the  horizontal 
four-bar  images  without  the  electronics  PTF  effects  at  1.15  cycles/mrad  including  and 
excluding  aliasing  effects  respectively.  The  net  result  of  aliasing  as  visually  observed  is 
that  the  two  bars  in  the  middle  seem  to  merge  into  a  wider  bar.  The  two  bars  are  no 
longer  present  at  their  expected  locations  in  the  image  and  the  four-bar  pattern  appears  as 
a  three-bar  pattern.  Figure  7.17  gives  the  plot  of  a  profde  passed  across  the  image  and 
shows  the  deterioration  of  the  symmetry  and  also  the  shift  in  the  locations  of  the  two  bars 
in  the  middle  due  to  aliasing  effects. 
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Figure  7.15.  Horizontal  Four-Bar  Image  without  Aliasing  and  without  Electronics 

PTF  at  1.15  Cycles/mrad. 
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Figure  7.16  Horizontal  Four-Bar  Image  with  Aliasing  and  without  Electronics 
PTF  at  1.15  Cycles/mrad.  Aliasing  Distorts  the  Symmetry  of  the  Image. 
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Figure  7.17.  A  Profile  across  the  Horizontal  Four-Bar  Pattern  without  Electronics 
PTF  at  1.15  Cycles/mrad.  Compare  with  Figures  7.15  and  7.16. 


The  following  Figures  7.18  and  7.19  illustrate  the  horizontal  four-bar  images  with 
the  electronics  PTF  effects  at  1.15  cycles/mrad  including  and  excluding  aliasing  effects 
respectively.  The  net  result  of  electronics  PTF  is  noticeable  as  degradation  in  contrast  on 
the  leftmost  bar.  Figure  7.20  gives  the  plot  of  a  profile  passed  across  the  image  and  shows 
the  contrast  degradation  on  the  leftmost  bar  and  the  deterioration  of  the  symmetry  in  the 
four-bar  pattern  due  to  electronics  PTF  effects.  The  entire  pattern  is  shifted  in  the  “down- 
scan”,  or  later  time  direction. 
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Figure  7.18.  Horizontal  Four-Bar  Image  without  Aliasing  and  with  Electronics 

PTF  at  1.15  Cycles/mrad. 
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Figure  7.19.  Horizontal  Four-Bar  Image  with  Aliasing  and  with  Electronics  PTF  at 
1.15  Cycles/mrad.  Compare  with  7.16.  PTF  Introduces  no  Apparent  Change  in  the 

Symmetry. 
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Figure  7.20.  A  Profile  across  the  Horizontal  Four-Bar  Pattern  with  Electronics 
PTF  at  1.15  Cycles/mrad.  Compare  with  Figure  7.17. 


E.  CONCLUSIONS 

The  results  obtained  from  the  experiments  with  the  Virtual  Thermal  Image- 
Processing  model  indicated  several  important  conclusions. 

•  Aliasing  effects  are  present  in  the  four-bar  pattern  images  at  spatial 
frequencies  below  the  Nyquist  limit.  Thus  modeling  aliasing  effects 
should  be  reconsidered  for  this  frequency  region  as  well  as  above  Nyquist. 

•  Aliasing  can  have  a  significant  visual  enhancing  effect  at  bar  spatial 
frequencies  less  than  the  Nyquist  limit 

•  Aliasing  effects  cannot  be  adequately  represented  as  noise  alone 

•  Thermal  imaging  devices  employing  different  sampling  frequencies  in 
horizontal  and  vertical  directions  experience  stronger  distortion  effects  due 
to  aliasing  in  the  direction  where  sampling  frequency  is  much  smaller 
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•  Aliasing  pattern  manifesting  itself  as  ripples  stays  steady  independent  of 
the  bar  spatial  frequency  in  all  output  images 

•  The  ripple  effect  introduced  by  aliasing  into  the  thennal  imagery  is  a 
function  of  spatial  sampling  frequency 

•  System  electronics  PTF  shows  its  effect  as  contrast  degradation  and 
deterioration  of  symmetry  that  results  in  the  shift  of  the  bar  locations 

This  work  and  its  conclusions  are  also  summarized  in  a  paper.  The  preprint  of  this 
paper,  which  is  accepted  for  presentation  at  the  35th  Asilomar  Conference  on  Signals, 
Systems,  and  Computers  on  November  4  -  November  7,  2001  is  provided  in  Appendix 
K. 

F.  DIRECTIONS  FOR  FURTHER  RESEARCH 

Experiments  with  the  Virtual  Thennal  Image-Processing  model  have  revealed 
several  advantages  of  its  operation  while  bringing  almost  endless  further  research 
opportunities.  From  these,  a  few  significant  topics  will  be  pointed  out  as  directions  for 
future  studies. 

First,  the  Virtual  Thermal  Image-Processing  model  can  be  improved  in  several 
ways.  The  current  version  produces  four-bar  thennal  images  without  taking  into  account 
the  MRTD  predictions.  The  model  should  be  applied  in  the  future  to  predict  subjective 
MRTD  from  the  virtual  images. 

Second,  the  model  currently  creates  thermal  images  of  only  horizontal  and 
vertical  bar  patterns.  In  the  future  it  should  be  able  to  produce  the  images  of  rotated  bar 
patterns  at  any  angles. 

Third,  the  Virtual  Thermal  Image-Processing  model  simulates  thermal  imaging 
systems  based  on  the  subjective  MRTD  measurement  scheme.  Ultimately,  the  observer 
decides  whether  the  bars  are  resolved  or  not  in  the  model  created  images.  However  the 
model  could  have  a  version  based  on  the  objective  MRTD  measurement  scheme  that 
leaves  the  observer  and  the  display  out  of  the  loop.  Therefore  the  model  could  be  used  to 
simulate  an  ATR  device. 
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Fourth,  the  model  could  include  a  more  complete  display  MTF  employing 
separate  MTFs  in  the  horizontal  and  vertical  dimensions  for  a  more  accurate  simulation 
of  the  thennal  imaging  systems. 

Fifth,  the  current  model  creates  the  input  target  pattern  only  once  and  then 
appropriately  scales  it  for  various  bar  spatial  frequencies.  However,  all  four-bar  patterns 
in  the  output  images  seem  to  have  the  same  physical  size  although  they  show  different 
scales.  To  correctly  present  the  bar  patterns  to  the  observer  with  respect  to  their  actual 
sizes,  the  model  could  create  an  individual  input  target  pattern  for  each  single  bar  spatial 
frequency. 

Finally,  the  Virtual  Thennal  Image-Processing  model  could  take  advantage  of  a 
more  user-friendly  interface  and  more  efficient  program  coding. 
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APPENDIX  A.  FLIR-92  SUB  SYSTEM  MTF  EQUATIONS  [REF.  2] 


A. 


PRE-FILTER  MTFS 

(a)  Optics  MTFs: 

(1)  Diffraction-limited  MTF 


HODL{fs)  = 


a  cos(  )  {[■ s ) 

1 

1 

1  ^ 

_ 1 

0.5 

D0  D„ 

L  Do  J 

(A.l) 


(2)  Geometric  Blur  MTF 

H0GB(fs)  =  e(-2^ 


(A.2) 


(b)  Detector  Spatial  MTF: 

sin  (nSJs) 

DsUs>  (nSzfs)  <A-3) 


(c)  Focal  Plane  Array  Integration  MTF: 

sinO? fTvsti) 


HMs)  = 


0?jW/) 


(d)  Sample-scene  Phase  MTF: 


Hssp(fs)  =  [cos 


\rNj 


0J 


(e)  Image  Motion  MTFs: 


(1)  Linear  Image  Motion  MTF 


(A.4) 


(A. 5) 
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B. 


sin^v^,) 

n  ML\J  s)  /  r  >  \ 

WsVJi) 

(2)  Random  Image  Motion  MTF 


HMR{fs)  =  e 


(A. 6) 


(A.  7) 


(3)  Sinusoidal  Image  Motion  MTF 


Hm<Js)  =  Jj?.*Afs) 


TEMPORAL  POST-FILTER  MTFS 

(a)  Detector  Temporal  MTFs: 


f  /  ,  \  2  ^/2 


1  + 


V 


rA' 

\/dt  J 


(b)  Electronics  Low  Frequency  Response: 


f  /  .  \2»V1/2 


HELP{fr)  = 


1  + 


/, 


J 


\  V  j  ELP  J 

(c)  Electronics  High  Frequency  Response: 


HEHP{fT)  = 


f  r  \n  (  f  - 


fr 


f 


/ 


V  J  EHP  J 


1  + 


fr 


^  V  f  EHP  J  J  j 


(d)  Boosting  MTF: 


HEB{fT)  =  1  + 


A-i) 


f 


1  -  COS 


V 


f 

V  fn  JJ 


(A. 8) 


(A. 9) 


(A.  10) 


(A.  11) 


(A.  12) 
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c. 


SPATIAL  POST-FILTER  MTFS 

(a)  Electro-optical  Multiplexor  MTF: 

sin  (xSledfs) 


H  MUX  if s  )  = 


Wledfs) 


(b)  Digital  F  ilter  MTF : 

(JV— 1)/2 

HDIG(fS)=  X  «,C0S 

1=0 


V  fco  J 


N/2 


f  1  ^ 


HDIG(fs)  =  'ZaiC0S 


i= 1 


fc 


CO 


(c)  CRT  Display  MTF: 


H CRT  (  f S  )  ^ 


<-2/rV/s2) 


a  = 


-log  0.025 

( N 

In1  r 

U; 


(A.  13) 


(A.  14) 


(A. 15) 


(A.  16) 


(A.  17) 


(1)  CCD  Charge  Transfer  Efficiency  MTF: 


Hccd{fs)  =  ^ 


-N(l-e )(l-cos  (2irfs  / rN  ))] 


(2)  Display  Sample  and  Hold  MTF: 


HDSH(fs) 


sin(  7idsfs ) 

(xSJs) 


(A.  18) 


(A.  19) 
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(A. 20) 


(d)  Eye  MTF: 

(1)  Non-limiting  Eye  MTF 

HEYE{fs)  =  1-0 

(2)  Limiting  Eye  MTF 

HEYE{fs)  =  e{2M]  (A-21) 
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APPENDIX  B.  DIRECTIONAL  AVERAGING  OPERATOR  [REF.  2] 


The  3-D  noise  concept  utilizes  an  averaging  process  to  isolate  the  noise 
components  in  the  desired  directions  using  the  directional  averaging  operators:  Dj,  Dr 
and  Dy.  Subscripts  of  these  operators  refer  to  the  directions  in  which  the  averaging 
operation  is  applied.  Each  D  operator  has  the  classical  simple  averaging  effect  in  the 
calculations.  The  original  three-dimensional  image  data  set,  which  contains  the  complex 
noise  components,  is  represented  by  U(T,V,H). 

The  three  important  properties  of  the  D  operators  are: 

a.  Commutative  Property: 

Dv  Dh  {U(T,V,H)}=  Dh  Dv  {U(T,V,H)} 

b.  Distributive  Property: 

Dv  {U(T,V,H)+  V(T,V,H)}=  Dv  {U(T,V,H)}+  Dv  {V(T,V,H)} 

c.  Idem-potent  Property: 

Dy  Dv  {U(T,V,H)}=  Dv  {U(T,V,H)} 

In  order  to  extract  the  desired  noise  component  from  the  composite  noise  data  set 
U  (T,  V,  H),  the  appropriate  directional  averaging  operators  should  be  used  as  shown  in 
the  table  below. 


NOISE  TERM 

CORRESPONDING  DIRECTIONAL 
OPERATION  REQUIRED 

Ntvh 

[(1-Dt)(1-Dv)(1-Dh)1  {U(T,V,H)| 

Nvh 

r (Dt)(  1  -Dy)(  1  -Dh)]  {U(T,V,H)| 

Ntv 

[(1-Dt)(1-Dv)(Dh)]  {U(T,V,H)} 

Nv 

[(DT)(l-Dy)(DH)l  |U(T,V,H)} 

Nth 

[(1-Dt)(Dv)(1-Dh)1  |U(T,V,H)} 

Nh 

[(Dt)(Dv)(1-Dh)1  |U(T,V,H)} 

Nt 

[(1-Dt)(Dv)(Dh)1  |U(T,V,H)} 

S 

r(DT)(Dv)(DH)l  {U(T,V,H)| 

Table  B.l.  The  Use  of  D  Operators  [Ref.  2], 
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The  D  operator  has  the  effect  of  deleting  the  selected  data  types  from  the  original 
composite  set.  On  the  other  hand  the  operations  (1-DT),  (1-DV),  and  (1-DH)  have  the 
opposite  effect.  They  extract  the  desired  noise  components  from  the  original  composite 
data.  For  example,  the  operation  [(Dt)(1-Dv)(Dh)]{U(T,V,H)}  indicates  that  the  desired 
component  to  be  extracted  from  the  composite  data  is  the  vertical  noise  component 
because  the  Dt  and  Dr  operators  remove  the  components  in  the  horizontal  and  temporal 
directions  while  the  (1-Dy)  operation  extracts  the  desired  vertical  noise  component.  The 
applications  of  each  directional  operator  are  also  presented  in  the  following  table. 


APPLIED  OPERATION 

DELETES 

EXTRACTS 

(Dt) 

Nj  ,  Njvh  ,  Njh  ,  Ntv 

S  ,  Nh  ,  Nvh  ,  Nv 

(1-Dt) 

S  ,  Nh  ,  Nvh  ,  Ny 

Nt  ,  Ntvh  ,  Nth  ,  Ntv 

(Dh) 

Nth  ,  Nr  ,  Nvh  ,  Ntvh 

S  ,  Nt  ,  Nv ,  Ntv 

(1-Dh) 

S  ,  Nt  ,  Ny ,  Ntv 

Nth  ,  Nh  ,  Nvh  ,  Ntvh 

(Dv) 

Ntvh  ,  Ny  ,  Nhv  ,  Ntv 

S  ,  Nh  ,  Nt,  Nth 

(1-Dv) 

S  ,  Nh  ,  Nx,  NTh 

Ntvh  ,  Ny  ,  Nhv  ,  Ntv 

Table  B.2.  The  Applications  of  Directional  Averaging  Operations. 
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APPENDIX  C.  EYE-BRAIN  SPATIAL  INTEGRATION  FACTORS 
FOR  MRTD  PREDICTIONS  [REF.  2] 


A.  SCANNING  SYSTEMS 

■  The  horizontal  eye-brain  integration  function  in  the  horizontal  direction  is  calculated 


A/p  o 


(C.l) 


The  vertical  eye-brain  integration  function  in  the  horizontal  direction  is  calculated  as: 


faivtv)  dv 


(C.2) 


The  horizontal  eye-brain  integration  function  in  the  vertical  direction  is  calculated  as: 


£»,  (/s)  =  ^T  JS(v)ff^(v) 

A/p  o 


(C.3) 


The  vertical  eye-brain  integration  function  in  the  vertical  direction  is  calculated  as: 


EVv{fs)  =  ~  )H2NFviy)  — ^ 

S  /TV 

^  i/  —rr 


(C.4) 
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■  The  parameters  used  in  Equation  (C.l)  through  Equation  (C.4)  are  described  in  the 
table  below. 


PARAMETER 

DESCRIPTION 

UNITS 

Vs 

Scan  velocity 

mrad/s 

Afp 

System  noise  bandwidth 

Hz 

S(v) 

Detector  noise  power 
spectrum 

— 

6V 

Detector  vertical 
instantaneous  field  of  view 
(IFOV) 

mrad 

Sv 

Samples  per  detector 
vertical  IFOV 

— 

Hnfh(v) 

Horizontal  system  noise 
filter  MTF 

H\iv(v) 

Vertical  system  noise  filter 
MTF 

Table  C.l.  Eye-Brain  Integration  Factor  Parameters  for  Scanning  Systems. 


B.  STARING  SYSTEMS 

■  The  horizontal  eye-brain  integration  function  in  the  horizontal  direction  is  calculated 


as: 


8  1  “ 

E„Jfs)  =  iL 


sin 


^  7TV  ^ 


7TV 

2 7s 


dv 


(C.5) 


The  vertical  eye-brain  integration  function  in  the  horizontal  direction  is  calculated  as: 


£,„(/s)=— 

s„ 


sin 


^  Itw^ 


-l2 


Itw 

2 Js 


dv 


(C-6) 
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The  horizontal  eye-brain  integration  function  in  the  vertical  direction  is  calculated  as: 


Su 


sin 


r  Itzv' 

Ws) 


Ittv 

2 7s 


dv 


(C.7) 


The  vertical  eye-brain  integration  function  in  the  vertical  direction  is  calculated  as: 


EAfs) 


S, 


J  Hi,  M 


sin 


2/, 


7TV 


2  fs 


dv 


(C.8) 


The  parameters  used  in  Equation  (C.5)  through  Equation  (C.8)  are  described  in  the 
table  below. 


PARAMETER 

DESCRIPTION 

UNITS 

SH 

Detector  horizontal 
instantaneous  field  of  view 
(IFOV) 

mrad/sec 

Sh 

Samples  per  detector 
horizontal  IFOV 

— 

5V 

Detector  vertical 
instantaneous  field  of  view 
(IFOV) 

mrad 

Sv 

Samples  per  detector 
vertical  IFOV 

— 

Hnfh(v) 

Horizontal  system  noise 
filter  MTF 

— 

Hnfv(v) 

Vertical  system  noise  filter 
MTF 

— 

Table  C.2.  Eye-Brain  Integration  Factor  Parameters  for  Staring  Systems. 
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APPENDIX  D.  THE  MITSUBISHI  ELECTRONICS  IR-M500 
THERMAL  IMAGING  SYSTEM  PARAMETERS 


In  the  subjective  and  objective  MRTD  laboratory  measurements  only  the 
Mitsubishi  Electronics  IR-M500  thermal  imager  was  used.  Therefore  the  same  thermal 
imaging  system  was  simulated  in  the  Visibility,  the  FLIR-92  and  the  Virtual  Thermal 
Image-Processing  models  utilizing  the  system  parameters  as  given  below  in  Table  D.l. 


Number  of  horizontal  detectors 

512 

Number  of  vertical  detectors 

512 

Detector  horizontal  dimension 

26  pm 

Detector  vertical  dimension 

20  pm 

Detector  horizontal  active  dimension 

16.24  pm 

Detector  vertical  active  dimension 

12.49  pm 

D* 

5x  1 0 1  °cmHz°'5/watt 

Frame  rate 

60Hz 

Spectral  cut-on 

3  pm 

Spectral  cut-off 

5  pm 

F/number 

1.4 

Focal  length  (optics) 

5  cm 

Transmittance  of  optics 

0.95 

Active  CRT  lines 

480 

3-D  noise  level 

Moderate 

Display  brightness 

lOmF 

Schottky  barrier  height 

22eV 

PtSi  emission  coefficient 

0.16  cV1 

Integration  time 

16145.833 

BFIP  performance 

Yes 

Table  D.l.  Mitsubishi  Electronics  IR-M500  Thermal  Imager  System  Parameters 

[Refs  .19  and  24], 
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APPENDIX  E.  MRTD  MEASUREMENT  RESULTS 


A.  SUBJECTIVE  MRTD  MEASUREMENT  RESULTS 


Spatial  frequency 
(cycles/mrad) 

Heat-up  MRTD 
results  (degrees  C) 

Cool-down 
MRTD  results 
(degrees  C) 

Average  MRTD 
results  (degrees  C) 

0.06 

0.11 

0.06 

0.08 

0.09 

0.22 

0.06 

0.14 

0.12 

0.28 

0.11 

0.19 

0.15 

0.28 

0.17 

0.22 

0.18 

0.28 

0.17 

0.22 

0.24 

0.33 

0.22 

0.28 

0.29 

0.39 

0.28 

0.33 

0.35 

0.44 

0.33 

0.39 

0.46 

0.50 

0.33 

0.42 

0.56 

0.56 

0.39 

0.47 

0.72 

0.72 

0.50 

0.61 

0.79 

0.83 

0.72 

0.78 

0.88 

1.06 

0.83 

0.94 

1.13 

2.06 

1.00 

1.53 

1.17 

2.83 

2.83 

2.83 

1.21 

8.00 

6.78 

7.39 

Table  E.l.  Horizontal  MRTD  Measurement  Results. 
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Spatial  frequency 
(cycles/mrad) 

Heat-up  MRTD 
results  (degrees  C) 

Cool-down 
MRTD  results 
(degrees  C) 

Average  MRTD 
results  (degrees  C) 

0.06 

0.21 

0.09 

0.15 

0.09 

0.29 

0.21 

0.25 

0.12 

0.29 

0.21 

0.25 

0.15 

0.40 

0.30 

0.35 

0.18 

0.42 

0.38 

0.40 

0.24 

0.61 

0.50 

0.55 

0.29 

0.68 

0.53 

0.65 

0.35 

0.74 

0.66 

0.70 

0.46 

0.76 

0.75 

0.75 

0.56 

1.40 

1.10 

1.25 

0.72 

4.25 

3.05 

3.65 

0.79 

6.20 

4.60 

5.40 

0.88 

10.40 

7.40 

8.90 

1.13 

24.00 

16.00 

20.00 

Table  E.2.  Vertical  MRTD  Measurement  Results. 
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Figure  E.l.  Horizontal  MRTD  Measurement  Results. 
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Figure  E.2.  Vertical  MRTD  Measurement  Results. 
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Spatial  Frequency 
(cycles/mrad) 

Horizontal  MRTD 
(degrees  C) 

Vertical 
MRTD 
(degrees  C) 

2-D  MRTD 
(degrees  C) 

0.06 

0.08 

0.15 

0.1095 

0.09 

0.14 

0.25 

0.1871 

0.12 

0.19 

0.25 

0.2179 

0.15 

0.22 

0.35 

0.2775 

0.18 

0.22 

0.40 

0.2966 

0.24 

0.28 

0.55 

0.3924 

0.29 

0.33 

0.65 

0.4631 

0.35 

0.39 

0.70 

0.5225 

0.46 

0.42 

0.75 

0.5612 

0.56 

0.47 

1.25 

0.7665 

0.72 

0.61 

3.65 

1.4921 

0.79 

0.78 

5.40 

2.0523 

0.88 

0.94 

8.90 

2.8924 

1.13 

1.53 

20.00 

5.5317 

Table  E.3.  Average  of  the  MRTD  Measurement  Results. 
B.  OBJECTIVE  MRTD  MEASUREMENT  RESULTS 


Spatial  frequency 
(cycles/mrad) 

MRTD  (degrees  C) 

0.06 

1.33 

0.12 

1.50 

0.24 

1.72 

0.35 

2.06 

0.46 

2.39 

0.56 

2.83 

Table  E.4.  MRTD  Results  for  SNR  of  6.0  (Constant  SNR  Case). 


SNR 

Constant  spatial  freq 

|uency  (cycles/mrad) 

0.12 

0.24 

0.46 

0.56 

2.0 

0.17 

0.44 

0.78 

1.00 

3.0 

0.39 

0.72 

1.33 

1.44 

4.0 

0.94 

1.11 

1.72 

1.83 

5.0 

1.22 

1.39 

2.17 

2.28 

6.0 

1.56 

1.78 

2.44 

2.89 

Table  E.5.  MRTD  Results  for  Constant  Spatial  Frequency  Case. 
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APPENDIX  F.  VISIBILITY  AND  FLIR-92  MODEL  MRTD 

PREDICTIONS 


Figure  F.l.  Visibility  Model  Subjective  MRTD  Predictions. 
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Figure  F.2.  FLIR-92  Model  MRTD  Predictions. 
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Figure  F.3.  Visibility  Model  Objective  MRTD  Predictions  for  SNR  of  6.0 

(Constant  SNR  Case). 


127 


Figure  F.4.  Visibility  Model  Objective  Predictions  MRTD  vs.  SNR  (Constant 

Spatial  Frequency  Case). 
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APPENDIX  G.  THE  VIRTUAL  THERMAL  IMAGE-PROCESSING 

MODEL  CODE 


The  original  MATLAB  code  files  that  are  used  in  generating  the  virtual  thermal 
images  are  contained  in  this  appendix.  “Virtual”  is  the  core  program  that  contains  the 
input  parameters  of  the  thermal  imager  system  being  modeled  and  calculates  the  total 
system  noise  variance.  Then  it  calls  a  second  program  “Vitualmtf’  to  calculate  the 
subsystem  MTFs  and  to  introduce  the  white  noise  into  he  system.  At  this  point,  the  white 
noise  is  generated  according  to  the  noise  variance  value  calculated  in  “Virtual”. 
“Vitualmtf’  then  calls  a  third  program  “alias”  to  include  the  aliasing  effects  into  the 
calculations. 


A.  VIRTUAL 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%  THE  VIRTUAL  THERMAL  IMAGE-PROCESSING  MODEL 


%  Yucel  Renter 

%  7/25/2001 

%  Given  system  parameters  (Mitsubishi  for  this  thesis), 

%  this  Matlab  file  will  calculate  the  total  system  noise 
%  variance  that  will  be  used  to  generate  the  white  noise 
%  Required  m-files: 

%  Virtual. m,  Virtualmtf.m,  alias. m,  Qs2var.m, 

%  Quad2var.m,  Sre.m 

%%%%%%0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 

clear; 

elf; 

0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 


%%%  Spatial  frequencies  of  interest  (cy/mrad)%%%%%%%%%%%%%% 
fbar  =  [.05:. 05:1.4]; 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
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%%%%%%%%%%%%%%%%%  System  Parameters  %%%%%%%%%%%%%%%% 
nh=5 12;  %  number  of  horizontal  detector  elements 


nv=512; 

%  number  of  vertical  detector  elements 

a=16.24; 

%  detector  active  horizonal  dimension  (um) 

b=12.49; 

%  detector  active  vertical  dimension  (um) 

hpitch=26; 

vpitch=20; 

%  detector  pitch,  horizontal  (um) 

%  detector  pitch,  vertical  (um) 

fnumber=  1.4;  %  system  f-number 
focal=50;  %  effective  focal  length  (mm) 

Fdot=60;  %  frame  rate  (Hz) 

clock=.97;  %  clock-out  factor  (for  integration  time) 
Dstarp=5el0;  %  peak  D-star  (cm -Hz'"  1/2/ Watt) 
lambdap=5.0;  %  peak  wavelength  (um) 
lambdal=3.0;  %  wavelength  band  lower  limit  (um) 
lambda2=5.0;  %  wavelength  band  upper  limit  (um) 


fxec=2.3e6; 

%  electronic  cut  off  (Hz) 

monfac=.25; 

%  monitor  Gaussian  rms  factor  normalized  by  alpha 

poles=l; 

blurspot=0; 

theh=0; 

%  number  of  poles  for  electronic  filter 

%  standard  deviation  of  blur  spot  diameter 

%  sample  scene  phase  angle,  horizontal  (rad) 

thev=0; 

%  sample  scene  phase  angle,  vertical  (rad) 

lloydc2  =  1.4388e4; 

%  c2  constant  in  Planck's  blackbody  equation  (used  in  delta  Tsc  calculation) 
Tbg=300;  %  background  temperature  (Kelvin) 


To=95; 

%  optical  transmission 

delTsc=0; 

%  minimum  threshold  input  contrast 

sigmavh=.4; 

SNRthr=6.0; 

%  fixed  pattern  noise  (percent  of  sigma  tvh) 

%  threshold  signal-to-noise  ratio 

te=.2; 

%  eye  integration  time 
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fpeye  =0.4;  %  eye  peak  response  frequency 

Deye  =3.3;  %  diameter  of  the  eye  pupil  (mm) 

monwave=0 .55;  %  monitor  wavelength  (microns) 

%0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 


%%%%%%%%%%%%%%%%%  Preliminary  Calculations  %%%%%%%%%%%%% 

%0/o0/o0/o0/o0/o%0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%0/o0/o0/o0/o0/o0/o%0/o0/o0/o0/o0/o%0/o0/o0/o%%% 


avglambda=0 . 5  *  (lambda  1 +lambda2); 

opdia=focal/fnumber; 

alpha=a/focal; 

beta=b/focal; 

sih=hpitch/focal; 

siv=vpitch/focal; 

ao=pi*(opdia/2)A2; 

td=clock/Fdot; 

fopt=opdia/avglambda; 

fsh=focal/hpitch; 

fsv=focal/vpitch; 

tclock=nh*nv*Fdot; 


fehis=fxec/fth; 


sspmaxh=0 
sspmaxv=0; 
feyec=Deye  /monwave; 


%  average  wavelength 
%  optical  diameter(mm) 

%  horizontal  DAS  in  mrad 
%  vertical  DAS  in  mrad 
%  horizontal  sampling  interval  in  mrad 
%  vertical  sampling  interval  in  mrad 
%  area  of  collecting  lens  (mmA  2) 

%  detector  dwell  time  (sec) 

%  optical  cut  off  in  cy/mr 
%  horizontal  sampling  frequency  (cy/mr) 

%  vertical  sampling  frequency  (cy/mr) 

%  in-clock  direction,  clock-out  frequency  (in  Hz) 
fth=sih*tclock;  %  clock- 

out  frequency  conversion  factor,  in-clock  direction 
ftv=siv*tclock/nh;  %  clock- 

out  frequency  conversion  factor,  out-of-  clock 
direction 

%  in-clock  direction,  electronic  cutoff  in  (cy/mrad) 
fevos=fxec/ftv;  %  out-of¬ 

clock  direction,  electronic  cutoff  in  (cy/mrad) 

%  sample  scene  phase  error,  horizontal 
%  sample  scene  phase  error,  horizontal 
%  eye  cut  off  in  cy/mrad 
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%%%0/o0/o0/o%0/o0/o0/o0/o0/o0/o0/o0/o0/o%0/o0/o0/o0/o0/o0/o%0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%0/o0/o%%%% 


%%%%%%%%%%%%%%%%%  Total  system  noise  calculations  %%%%%%%%%% 

%%%%%%0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 

[MTFrf,fmax]= 

meffcnt(fbar,fxec,fsh,fsv,sih,siv,fth,ftv,fopt, alpha, beta, blurspot,theh,thev,fehis,fevos, pole 
s,monfac,sspmaxh,sspmaxv,feyec,fpeye); 


%%%%%%%%%%  NETD  calculation  %%%%%%%%%%%%%%%%%%%%%%% 

refbw  =  (pi/4)*(  1/td);  %  reference  bandwidth  calculation 

nl  =  pi*sqrt(a*le-4*b*le-4*refbw); 

dl  =  (alpha*beta*le-6*To*ao*le-2*Dstarp); 

d2  =  lloydc2/(lambdap*TbgA2)*(quad2var('sre',  lambda  1,  lambda2,  Tbg)); 

NETD  =  nl/(dl*d2) 

%%%%%%%%%%  SNR  improvement  calculation  %%%%%%%%%%%%%%%%%% 

tn2  =  sqrt(  1 +sigmavh  A2/ (te*Fdot))*(piA2*  SNRthr/ ( 8  *  sqrt(te  *  Fdot))) 

%%%%%%%%%%  Total  system  noise  calculation  %%%%%%%%%%%%%%%%% 

MTFsys  =  (abs(MTFrf)).A2;  %  overall  system  MTF 

totbw  =  sum(sum(MTFsys)).*(alpha/td).*(fmax/128) 

bwcorr  =  sqrt(totbw/refbw)  %  bandwidth  correction  factor 

sigmatvh=NETD*bwcorr; 

sigmatotal=(sigmatvhA2+(0 .4  *  sigmatvh)A2)A .  5 ; 

vartotal  =  sigmatotalA2  %  Total  noise  variance  (as  in  FLIR-92) 

%%%%  Correction  to  noise  variance  for  compensation  %%%%%%%%%%%%%%% 
varcorrected  =  vartotal*(fmaxA2/totbw) 

%Same  value  is  used  in  meffcnt.m  for  the  generation  of  white  noise 

0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 


B.  VIRTU  ALMTF 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%  Yucel  Renter 


%  This  program  creates  virtual  standard  four-bar  patterns 
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%  as  input  to  the  simulated  TIS  (Mitsubishi). 

%  It  then  calculates  the  subsystem  MTFs  and  applies  them  to 
%  the  input  and  finally  produces  the  four-bar  pattern 
%  images  at  the  output  for  the  spatial  frequencies 
%  specified  by  the  user.  Overall  system  (Gaussian  white) 

%  noise  is  also  created  and  added  here  at  the  detector. 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


function 


[MTFrf, fmax]=meffcnt(fbar,fxec,fsh,fsv,sih,siv,fth,ftv,fopt, alpha, beta, blurspot,theh,thev, 

fehis,fevos,  poles,  monfac,sspmaxh,sspmaxv,feyec,fpeye); 

d2=length(fbar); 

0/o0/o0/o0/o%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%%%%  Create  bar  pattern  (Input  contrasts)  %%%%%%%%%%%%%%% 
coll  =  zeros(130,10); 

col2  =  [zeros(30,10);ones(70,10);zeros(30,10)]; 
col3  =  zeros(  130,93); 
rowl  =  zeros(63,256); 

%Horizontal  MRTD  bar  pattern  (A) 

A=[rowl  ;col3  ,col2,col  1  ,col2,col  1  ,col2,col  1  ,col2,col3  ;row  1  ] ; 

%  Vertical  MRTD  bar  pattern  (V) 

V=A'; 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%%%%%  Plot  bar  pattern  at  the  input  %%%%%%%%%%%%%%%%% 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


xside=linspaee(0, 12,73); 

yside=xside; 

%HORIZONTAL 


figure(l) 

surf(xside,yside,abs(A(56:2:200, 56:2:200))) 
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shading  interp; 
axis([0  12  0  12  0  2]) 
axis  off 

xlabel('x  space,  mrad') 
ylabel('y  space,  mrad’) 
colonnap(gray) 
view(2) 

%VERTICAL 


figure(2) 

surf(xside,yside,abs(V(56:2:200, 56:2:200))) 
shading  interp; 
axis([0  12  0  12  0  2]) 
axis  off 

xlabel('x  space,  mrad’) 
ylabel('y  space,  mrad') 
colormap(gray) 
view(2) 

%0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 


%%%%%%%%%%  Transform  into  spatial  freq.  domain  %%%%%%%%%%%%%%% 
%%%%%%%%  and  apply  scaling  to  represent  the  appropriate  %%%%%%%%%%%% 
%%%%%%%%%%%%%%  spatial  frequency  %%%%%%%%%%%%%%%%%%%% 
B=  fft2(A);w=fft2(V); 
b=abs(B);x=abs(w); 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%%%%%%%%%%%  Create  Gaussian  white  noise  %%%%%%%%%%% 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


vartotal  =  0.0022;  %  This  value  is  the  calculated  total 

noise  value  for  delTsc=0.0981  (vartotal=sigmatvhA2  + 
sigmavhA2) 
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varcorrected  =  0.0087;  %  This  value  is  the  corrected  total  noise  value 

noise  =  sqrt(varcorrected).*randn(size(B));  %  Noise  matrix  in  space 

N  =  fftshift(1.45.*fft2(noise));  %  Noise 
spectrum  (1.45  is  the  compensation  factor  for  the 
contrast  manipulation  accomplished  after 
reconstruction  filtering) 

%%%%%%%%%%%%%%%%%  Now  loop  through  frequencies  %%%%%%%%%% 

for  counter  =  l:d2; 
fx  =  fbar(counter); 

0/o0/o%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%  W  refers  to  the  actual  width  of  one  bar  in  mrad  %%%%%%%%%%%%% 
%%%%%%%  'del'  terms  are  used  to  find  the  appropriate  %%%%%%%%%%%%%%% 
%%%%%%%%%%%  frequency  and  space  scales  %%%%%%%%%%%%%%%%%% 
W=l/(2*fx); 

Nw  =10;  %  number  of  elements  in  one  bar  (x) 

Nt  =  256;  %  number  of  elements  in  vector  (m  x  m) 

delt  =  W/Nw; 

DELt  =  Nt*W/Nw; 
delf  =  1/DELt; 

DELf  =  1/delt; 
fimax  =  DEL172; 

fscale  =  linspace(-fmax,  fimax, 256); 
xscale  =  linspace(-DELt/2,DELt/2,256); 
yscale  =  linspace(-DELt/2,DELt/2,256); 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%%%%  Create  a  spatial  freq.  matrix  to  be  used  in  %%%%%%%%%%%% 
%%%%%%%%%%%%%  MTF  calculations  %%%%%%%%%%%%%%%%%%%% 

uux=l:256; 

uuy=l:256; 


135 


Matfx=uux*uuy'; 
Matfx=ones(size((Matfx))); 
Matfy=Matfx; 
for  k=- 128  :  127; 
for  n=-128:127; 
k2=k+129; 


n2=n+129; 

Matfx(k2  ,n2)=k*  fmax/ 128; 
Matfy(k2,n2)=n*fmax/ 128; 
end; 
end; 

fr=((Matfx).A2+(Matfy).A2).A(0.5); 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%%%%%%%%%%%  Create  image  formation  MTF  %%%%%%%%%%% 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%  DIFFRACTION  LIMITED  OPTICAL  MTF. 


(Circular  aperture) 
cl=fr/fopt; 
for  ii=l:256; 
for  jj=l:256; 
if  cl(ii,jj)>1.0; 

cl(iijj)=1.0; 

end; 

end; 

end; 

c2=sqrt(l-cl  .A2); 

ro=2/pi  .*(acos(cl)-(cl  .*c2)); 

%  GEOMETRIC  BLUR  OPTICAL  MTF 
fl  =  (fr.A2).*piA2.*(-2).*blurspotA2; 
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rb  =  exp(fl); 


%  HORIZONTAL  AND  VERTICAL  SAMPLE  SCENE  PHASE  MTF 

rsh  =  (2*Matfx/fsh)*theh; 

rsh2  =  cos(rsh); 

rsv  =  (2*Matfy/fsv)*thev; 

rsv2  =  cos(rsv); 

rs  =  rsh2.*rsv2; 

%  DETECTOR  HORIZONTAL  MTF.  (For  rectangular  detector) 
c3fx=pi.*alpha.*Matfx; 

%  DETECTOR  VERTICAL  MTF.  (For  rectangular  detector) 
c3fy  =pi.*beta.*Matfy; 
for  m=l  :256; 
for  p=  1:256; 
if  c3fx(m,p)==0 
c3fx(m,p)=0.01; 
end; 

if  c3fy(m,p)==0 
c3fy(m,p)=0.01; 
end; 
end; 
end; 

rdh=abs(sin(c3fx)./c3fx)  ; 
rdv=abs(sin(c3fy)./c3fy) ; 

%  DETECTOR  2D  MTF 
rd=rdv.*rdh; 

%  Cascade  MTFs 
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u=ro.*rd.*rb.*rs; 


0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%%%%%%%%%%%%%%%%%%%%% 

%%%%%%  Add  white  noise  and  multiply  image  spectrum  by  image  %%%  %%%%% 
%%%%%%%%%%%%%%%%formation  MTF  %%%%%%%%%%%%%%%%%%% 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


d=fftshift(B);dv=  fftshift(w);display('showing  noise  term');fbar(counter) 
mess  =  d.*u+N;mes2  =  dv.*u+N;  %  noise  is  added  here 

0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 


%%%%%%%%%%%%%%%%%  Calculate  aliasing  terms  %%%%%%%%%%%%%% 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


aliasterm  =  zeros(size(mess)); 
valiasterm  =  aliasterm; 

[aliasterm,  aliasx]  = 

alias(Matfx,Matfy,fscale,fsh,fsv, mess, sspmaxh,sspmaxv,fbar, counter); 

[valiasterm,  aliasx]  = 

alias(Matfx,Matfy,fscale,fsh,fsv,mes2,sspmaxh,sspmaxv,fbar,  counter); 

0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 


%%%%%%%%%%%%  Add  alias  tenn  to  pre-filtered  scene  %%%%%%%%%%%%% 

0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 


withalias=  mess  +  aliasterm; 
vertalias  =  mes2  +  valiasterm; 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%%  Create  image  reconstruction  MTF  %%%%%%%%%%%%%%%%%% 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%  ELECTRONICS  OTF.  (Multi-pole  low  pass  filter) 
%  OTF  from  impulse  response  function) 

N=256; 

dx=l/(2*fmax); 
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dy=dx; 

uux2=(0:(N-l))*dx; 

uuy2=(0 :  (N- 1 ))  *  dy ; 

respx=exp(-uux2  *  fehis*  (2  *pi)); 

respy=exp(-uuy2  *  fehis*  (2  *pi)); 

respxy=respx'  *  respy; 

sgnfhp=((-  l)*ones(size(  1  :N))).A(0:(N- 1)); 

uuxy=sgnfhp'*sgnflip; 

respxy2=respxy .  *uuxy; 

eloTF=FFT2(respxy2); 

mag=abs(eloTF); 

eloTF=eloTF/(max(max(mag))); 

%  Nonnalization  is  done  here 
elMTF=abs(eloTF); 

%  CRT  MONITOR  SPATIAL  MTF. 
c7=(fr/(l/sih))  ,A2; 

c8=(monfac)A2;  %  sigma/sih  =  .25 
rm=exp(-2  *piA2  ,*c7  .*c8); 

%  CASCADE  ELECTRONICS  OTF  AND  DISPLAY  MTF 
recon=rm.  *eloTF.  *rd; 

%0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 


%%  Calculate  system  MTF  (for  noise  bandwidth  calculations) 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%  Calculate  eye  MTF 
cl  l=fr/feyec; 
for  iii=l:256; 


for  jjj= 1:256; 
if  cl  l(iii,jjj)>1.0; 
cl  l(iii,jjj)=1.0; 
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end; 

end; 

end; 

c22=sqrt(  1  -c  1 1 .  A2); 
constl=0.4/feyec; 
const2=sqrt(  1  -const  1 A2); 

nonn=  l/(2/pi  *(acos(constl)-(constl*const2))); 

c33=2*norm/pi  ,*(acos(cl l)-(cl 1  ,*c22)); 

clear  cl  1  c22; 

c44=fr./fpeye; 

reye=min(c44,c33); 

clear  ell  c22  c33  c44; 

%  Multiply  Eye  MTF  with  the  Reconstruction  OTF 
MTFrf  =  reye.*recon; 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%  Multiply  pre-fdtered  and  aliased  image  spectrum  %%%%%%%%%%%%%%% 
%%%%%%%  with  the  aliased,  fdtered,  scene  %%%%%%%%%%%%%%%%%%%% 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


outalias  =  withalias.*recon; 
noalias  =  mess.*recon; 
veral2  =  vertalias.*rccon; 
veraB  =  mes2.*recon; 


%  Output  spectrum  WITH  aliasing  (Horizontal) 

%  Output  spectrum  WITHOUT  aliasing  (Horizontal) 
%  Output  spectrum  WITH  aliasing  (Vertical) 

%  Output  spectrum  WITHOUT  aliasing  (Vertical) 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%%%  Transform  back  to  spatial  domain  %%%%%%%%%%%%%%%%% 
%%%%%%%%%  A  normalization  function  is  used  %%%%%%%%%%%%%%%%% 
%%%%%%%%%  for  a  better  visualization  %%%%%%%%%%%%%%%%%%% 
%%%%%%%%%  normfunc  =  ((2/pi).*atan(5.*(x-0.5))+l)./2;  %%%%%%%%%%%% 

%0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 


messy  =  fftshift(  outalias); 
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messy2  =  fftshift(noalias); 
gigo  =  fftshift(veral2); 
gigo2  =  fftshift(veraB); 
C=ifft2(messy); 

C=((2/pi).*atan(5.*(C-0.5))+l)./2; 

D=ifft2(messy2); 

D=((2/pi).*atan(5.*(D-0.5))+l)./2; 

E=ifft2(gigo); 

E=((2/pi) .  *  atan(5 .  *  (E-0.5))+ 1  )./2; 
F=ifft2(gigo2); 

F=((2/pi).  *  atan(5 .  *  (F-0 .5))+ 1)./2; 


%  Output  image  WITH  aliasing  (Horizontal) 

%  Normalization  is  applied  here 
%  Output  image  WITHOUT  aliasing  (Horizontal) 
%  Normalization  is  applied  here 
%  Output  image  WITH  aliasing  (Vertical) 

%  Normalization  is  applied  here 
%  Output  image  WITHOUT  aliasing  (Vertical) 

%  Nonnalization  is  applied  here 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%%%%%%%%%  Virtual  images  are  created  here  %%%%%%%%%%%% 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%HORIZONTAL 

figure(3) 

xside=linspace(0,(145/(fx*20)),73);yside=xside; 
%  Modified  for  the  scaling  on  the  axis 
surf(xside,yside,abs(C(56:2:200, 56:2:200))) 
shading  interp; 

axis([0  (145/(fx*20))  0  (145/(fx*20))  0  2]) 

%  Modified  for  the  scaling  on  the  axis 

xlabel('x  space,  mrad') 

ylabel('y  space,  mrad') 

colormap(gray) 

grid  off 

view(2) 

%VERTICAL 

figure(4) 
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xside=linspace(0,(145/(fx*20)),73); 

%  Modified  for  the  scaling  on  the  axis 
yside=xside; 

surf(xside,yside,abs(E(56:2:200, 56:2:200))) 
shading  interp; 

axis([0  (145/(fx*20))  0  (145/(fx*20))  0  2]) 

%  Modified  for  the  scaling  on  the  axis 

xlabel('x  space,  mrad') 

ylabel('y  space,  mrad') 

colonnap(gray) 

grid  off 

view(2) 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


C.  ALIAS 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%  This  program  calculates  the  alias  tenns  and  turns  them 
%  back  to  “virtualmtf  ’ 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

function  [aliasterm,aliasx]  =  alias(Matfx,Matfy,fscale,fsh,fsv, mess, b,c,fbar, counter) 

0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o0/o%%% 


%%%%%%%  Represent  aliasing  by  repeating  the  filtered  %%%%%%%%%%%%%% 
%%%%%%%  spectrum  at  sampling  fireq.  intervals  %%%%%%%%%%%%%%%%%% 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 


%%%%%%%%%%  Alias  terms  (other  than  cross  aliases)  %%%%%%%%% 
where  =  find(fscale/fsh  >1); 
if  isempty(where); 

cent  =  256;multiple=0; 
else 
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fbar(counter) 
cent  =  where(l) 
multiple=l; 
end 

start  =  cent-128;unused  =  256-start; 

%%%%  Alias  terms  due  to  sampling  in  x  direction  %%%%% 

%First  positive  side 

aliastxp  =  circshift(mess,[0  (start- 1)]); 

aliastxp  =  aliastxp.  *exp(j*2*pi*b*fsh); 

%Now  negative  side 

aliastxn  =  circshift(mess,[0  -(start-1)]); 

aliastxn  =  aliastxn.*exp(-j*2*pi*b*fsh); 

%Combine  the  two 
aliasx  =  aliastxp+aliastxn; 
wherev  =  find(fscale/fsv  >1); 
if  isempty(wherev); 
centv  =  256; 
multiplev  =  0; 

else 

fbar(counter) 
centv  =  wherev(l) 
multiplev  =  1 ; 
end 

startv  =  centv- 128; 
unusedv  =  256-startv; 

%%%%  Alias  terms  due  to  sampling  in  y  direction  %%%%% 

%First  positive  side 

aliastyp  =  circshift(mess, [(startv- 1)  0]); 

aliastyp  =  aliastyp.  *exp(j*2*pi*c*fsv); 
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%Now  negative  side 

aliastyn  =  circshift(mess,[-(startv-l)  0]); 

aliastyn  =  aliastyn.*exp(-j*2*pi*c*fsv); 

%Combine  the  two 
aliasy  =  aliastyp+aliastyn; 

%%%%%%%%%%  Cross  alias  terms  %%%%%%%%%%%%%%%%%%%%%%% 

aliasxpyp=circshift(mess,[(startv- 1 )  (start- 1 )] ); 
aliasxpyp=aliasxpyp .  *  exp(j  *2*pi*b*  fsh) .  *  exp(j  *2*pi*c*  fsv); 


aliasxpyn=circshift(mess,[-(startv-l)  (start-1)]); 
aliasxpyn=aliasxpyn.*exp(j*2*pi*b*fsh).*exp(-j*2*pi*c*fsv); 


aliasxnyp  =  circshift(mess,[(startv-l)  -(start- 1)]); 
aliasxnyp=aliasxnyp.*exp(-j*2*pi*b*fsh).*exp(j*2*pi*c*fsv); 


aliasxnyn=circ shift(mess, [-(startv- 1 )  -(start- 1 )]) ; 
aliasxnyn=aliasxnyn.*exp(-j*2*pi*b*fsh).*exp(-j*2*pi*c*fsv); 


crossterm=multiple*multiplev*(aliasxpyp+aliasxpyn+aliasxnyp+aliasxnyn); 
%%%%%%%%%%%%  Total  alias  term  %%%%%%%%%%%%%%%%%%%%%%% 
aliasterm  =  multiple*aliasx  +  multiplev*aliasy  +crosstenn; 


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
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APPENDIX  H.  THE  IMAGE  RECONSTRUCTION  DETECTOR 

MTF  ANALYSIS 

As  mentioned  in  Chapter  IV,  both  the  Visibility  model  and  the  Virtual  Thermal 
Image-Processing  model  include  a  detector  MTF  into  the  image  reconstruction  process  as 
treated  by  Lloyd  in  Reference  10.  In  this  appendix,  the  analysis  of  the  detector  MTF  in 
the  backend  filtering  process  will  be  made  in  the  images  generated  by  Virtual  Thennal 
Image-Processing  model. 

The  experiments  conducted  with  the  Virtual  Thermal  Image-Processing  model 
showed  that  some  of  the  distortion  appearing  as  patterns  of  ripples  due  to  aliasing  was 
smoothed  out  from  the  image  with  the  introduction  of  the  detector  MTF  into  the 
reconstruction  process.  Figures  H.l  and  H.2  demonstrate  the  patterns  of  distortion  on  the 
horizontal  and  the  vertical  four-bar  targets  at  1.15  cycles/mrad  respectively  when  the 
detector  MTF  is  excluded.  Following  Figures  H.3  and  H.4  show  how  the  detector  MTF 
filters  out  some  of  these  aliasing  effects  from  the  image  when  it  is  included  into  the 
reconstruction  MTF. 
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x  space,  mrad 


Figure  H.l.  Horizontal  Four-Bar  Pattern  Image  Generated  wthout  the  Additional 
Detector  MTF  at  1.15  Cycles/Mrad.  Note  Ripple  Detail  in  the  Vertical  Direction. 


_ 

0  1  2  3  4  5  6 

x  space,  mrad 

Figure  H.2.  Vertical  Four-Bar  Pattern  Image  Generated  without  the  Additional 

Detector  MTF  at  1.15  Cycles/Mrad. 
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x  space,  mrad 


Figure  H.3.  Horizontal  Four-Bar  Pattern  Image  Generated  with  the  Additional 
Detector  MTF  at  1.15  Cycles/Mrad.  Note  Suppression  of  Vertical  Ripple 

Modulation. 


x  space,  mrad 


Figure  H.4.  Vertical  Four-Bar  Pattern  Image  Generated  with  the  Additional 

Detector  MTF  at  1.15  Cycles/Mrad. 
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The  most  significant  observation  to  be  made  from  these  figures  is  the  visual 
difference  of  the  filtering  process  employed  by  the  image  reconstruction  detector  MTF  in 
the  output  images.  Although  the  same  MTF  is  used  for  both  the  horizontal  and  vertical 
bar  patterns  the  detector  MTF  filtering  effect  in  the  horizontal  bar  pattern  image  is  much 
more  apparent  because  the  filtered  ripples  were  brought  about  by  the  aliasing  in  the 
vertical  direction  that  was  previously  confirmed  to  have  weaker  effects.  This  explanation 
still  applies  for  the  vertical  bar  pattern  image,  however  the  net  effect  is  not  as  easily 
observable  because  the  vertical  bar  pattern  is  strongly  distorted  by  aliasing  in  the 
horizontal  direction  such  that  the  bars  appear  to  be  dissected. 


148 


APPENDIX  I.  RESULTS  OBTAINED  FROM  THE  ANALYSIS  OF 
RIPPLE  EFFECTS  BROUGHT  ABOUT  BY  ALIASING 


/s 

(cycles/mrad) 

n 

(n.  /s 
) 

0.2 

33 

6.6 

0.3 

21 

6.3 

0.4 

16 

6.4 

0.5 

13 

6.5 

0.6 

10 

6 

0.7 

9 

6.3 

0.8 

8 

6.4 

0.9 

7 

6.3 

1 

6 

6 

1.1 

5 

5.5 

1.2 

5 

6 

Table  1.1.  Ripple  Effect  Analysis  Results. 


Where 

/ s  is  the  bar  spatial  frequency  in  cycles/inrad 

n  is  the  number  of  ripples  observed  along  a  bar 


The  multiplication  of  number  of  ripples  with  the  bar  spatial  frequency  stays  at  an 
approximately  constant  value  as  seen  in  the  rightmost  column  in  Table  1.1.  This  is  an 
indication  of  the  fact  that  the  pattern  introduced  by  aliasing  into  the  thermal  images  is 
stationary  for  all  bar  spatial  frequencies.  An  example  illustrating  the  number  of  ripples 
along  a  bar  from  a  vertical  bar  pattern  for  a  spatial  frequency  of  0.90  cycles/mrad  is  given 
below  in  Figure  1.1. 
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Figure  1.1.  An  Example  Illustrating  the  Ripple  Effect  Along  a  Vertical  Bar  at 
0.90  Cycles/Mrad.  This  Pattern  Remains  Constant  as  the  Spatial  Frequency  of  the 

Bar  Pattern  Changes. 
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APPENDIX  J.  THE  RESULTS  OBTAINED  FROM  THE  ANALYSIS 
FOR  DERIVING  SAMPLING  FRQUENCY  FROM  THE  RIPPLE 
EFFECTS  BROUGHT  BY  ALIASING 


/s 

(cycles/mrad) 

W=l/(2./s) 

L=7.W 

n 

Ripple  frequency 
=(n/L) 

0.2 

2.5 

17.5 

33 

1.885 

0.3 

1.667 

11.667 

21 

1.8 

0.4 

1.25 

8.75 

16 

1.828 

0.5 

1 

7 

13 

1.857 

0.6 

0.833 

5.833 

10 

1.714 

0.7 

0.714 

5 

9 

1.8 

0.8 

0.625 

4.375 

8 

1.828 

0.9 

0.556 

3.889 

7 

1.8 

1 

0.5 

3.5 

6 

1.714 

1.1 

0.454 

3.181 

5 

1.571 

1.2 

0.4167 

2.9167 

6.6 

1.714 

Sampling  frequency  (calculated 

1  from  system  parameters)  =  1.92  cycles/mrad 

Table  J.l.  Derivation  of  Sampling  Frequency  from  Ripple  Effects  Brought  by 

Aliasing. 


Where 

/ s  is  the  bar  spatial  frequency  in  cycles/mrad 

W  is  the  angular  subtense  on  a  single  bar  in  mrad 

L  is  the  length  of  one  bar  in  mrad 

n  is  the  number  of  ripples  observed  along  a  bar  at  the  corresponding 

bar  spatial  frequency 

(n/L)  gives  the  ripple  frequency  in  ripples/mrad 
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As  seen  in  the  rightmost  column  in  Table  J.l,  the  ripple  frequencies  derived  for 
each  bar  spatial  frequency  from  the  ripple  pattern  produced  by  aliasing  agrees  within  a 
reasonable  range  of  error  with  the  sampling  frequency  calculated  from  the  Mitsubishi 
system  parameters. 
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Abstract 

A  new  virtual  thermal  image-processing  model  is  introduced  in  this  paper.  This 
visualization  program  is  based  on  earlier  modeling  work  focused  on  predicting  the 
minimum  resolvable  temperature  (MRTD),  which  is  a  standard  performance  measure  for 
forward  looking  infrared  radar  (FLIR)  imaging  systems.  Relevant  filtering,  noise  and 
sampling  processes  are  included  in  the  visualization  model.  In  this  paper  it  will  be 
demonstrated  and  explained  that  aliasing  effects  in  thermal  images  of  four-bar  patterns 
cannot  in  general  be  adequately  modeled  as  noise.  In  particular,  the  simulation 
experiments  demonstrate  that  aliasing  can  have  a  noticeable  visual  enhancing  effect  at 
spatial  bar  frequencies  less  than  the  Nyquist  limit. 

1.  Introduction 

Minimum  resolvable  temperature  difference  (MRTD)  is  a  standard  performance 
measure  for  FLIR  thermal  imaging  systems. [Ref  1]  It  is  defined  as  the  temperature 
difference  between  a  four-bar  target  and  its  uniform  background,  which  is  required  by  a 
trained  observer  to  just  resolve  all  four  bars.  MRTD  is  a  function  of  four-bar  target  spatial 
frequency  as  demonstrated  in  Figure  K.l.  Spatial  frequency  is  defined  as: 


fo  = 


R 

2  W 


cy  /  mrad 


(1) 


where  W  is  bar  width  in  millimeters  and  R  is  target  range  in  meters. 
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2.  Sampling  and  Aliasing 

Current  state-of-the  art  thermal  imaging  systems  incorporate  staring  focal  plane 
arrays  (FPA)  which  sample  the  scene  spatially.  Figure  K.2  demonstrates  two-dimensional 
spatial  sampling  by  detector  elements  in  a  staring  FPA.  Horizontal  spatial  sampling 
frequency  is  calculated  as: 

f  cyhnrad  (2) 

Ax 

where  foptiCs  is  optical  focal  length  in  mm  and  Ax  is  pixel  pitch  in  pm. 


Sampling  and  aliasing  effects  are  more  easily  appreciated  when  examined  in  the 
spatial  frequency  domain.  Multiplication  in  the  space  domain  corresponds  to  convolution 
in  the  frequency  domain.  Therefore,  the  sampling  effect  produces  replication  of  the  object 
spectrum  at  integer  multiples  of  the  sampling  frequency. 


Sampling  theory  suggests  that  a  frequency  component  in  the  original  image  above 
the  Nyquist  limit  appears  as  a  lower  frequency  component  after  sampling  due  to  the 
effect  called  ‘aliasing’.  Equation  3  defines  the  Nyquist  frequency  limit.  The  effects  of 
aliasing  and  the  degree  to  which  it  interferes  with  recognition  depend  on  the  scene  that  is 
viewed.  Aliasing  tends  to  be  apparent  as  image  distortion  and  degradation  when  viewing 
periodic  targets  such  as  a  four  bar  pattern. 

fN  =  —  cyhnrad  (3) 

Different  models  have  been  proposed  to  predict  the  perfonnance  of  staring 
thennal  imaging  systems.  The  infrared  community  standard  FLIR92  model  predicts 
MRTD  results  below  the  Nyquist  limit  and  ignores  aliasing  effects.  A  visibility  model, 
proposed  by  R.J.Pieper  and  A.W.  Cooper  of  the  Naval  Postgraduate  School  in  1994,  also 
predicts  the  MRTD  of  staring  thermal  imaging  systems.  The  model  is  based  on  a 
minimum  threshold  input  contrast  parameter  and  a  contrast  reduction  factor  due  to 
aliasing  and  blurring  effects  [Ref.  2],  It  provides  MRTD  predictions  beyond  the  Nyquist 
frequency  limit.  MRTD  predictions  from  the  visibility  model  have  been  shown  to  be  in 
better  agreement  with  the  laboratory  measurements  [Ref.  3], 

3.  Virtual  thermal  image  processing  model 

Virtual  MRTD  experiments  have  been  perfonned  at  the  U.S.  Army  Night  Vision 
Electronics  Sensors  Directorate  using  simulations  based  on  the  FLIR92  model.  In  the 
experiments,  effects  of  blurring  and  noise  on  the  image  have  been  reported  [Ref.  4], 

A  new  Virtual  thennal  image  processing  model  (or  virtual  thermal  image 
processing  model)  is  introduced  in  this  paper.  The  model  is  based  on  the  Visibility 
MRTD  model  and  has  been  developed  at  the  Naval  Postgraduate  School.  It  takes  the 


154 


current  thermal  imaging  system  modeling  concerns  such  as  sampling  and  aliasing  into 
account.  From  the  basic  staring  imaging  system  parameters,  the  model  creates  the  system 
response  and  applies  it  to  the  input  four-bar  target.  The  model  provides  visual  images  that 
can  be  obtained  with  the  actual  imaging  system  being  modeled.  This  allows  the  user  to 
virtually  evaluate  the  effects  of  noise,  sampling  and  aliasing  on  the  imagery.  The  virtual 
thermal  image  processing  model  is  written  using  MATLAB  computational  software.  The 
model  block  diagram  is  presented  in  Figure  K.3. 


First,  a  two-dimensional  four-bar  target  is  generated.  The  model  works  with  both 
horizontal  and  vertical  bar  patterns.  Taking  the  two-dimensional  Fast  Fourier  Transform 
(2-D  FFT)  of  the  bar  pattern,  a  four-bar  target  spatial  frequency  spectrum  is  obtained. 
The  target  spectrum  is  then  multiplied  by  the  Image  Formation  Modulation  Transfer 
Function  (MTF).  At  this  stage,  white  noise  is  generated  and  added  to  the  spectrum. 
Following  the  application  of  sampling  and  aliasing  effects,  the  filtered  and  aliased 
spectrum  is  multiplied  by  the  image  reconstruction  MTF.  Finally,  taking  the  inverse  2-D 
FFT  gives  the  visual  four-bar  target  representation. 

4.  Aliasing  Effects  Below  the  Nyquist  Limit 

Experiments  with  the  virtual  thermal  image-processing  model  have  provided 
some  useful  insights  into  the  problems  encountered  in  modeling  aliasing  effects.  Aliasing 
effects  are  shown  to  be  present  in  the  image  of  four-bar  targets  at  spatial  frequencies 
below  the  Nyquist  limit.  The  spatial  sampling  frequency  from  the  system  parameters  is 
1.92  cycles/mrad.  Profiles  across  a  four-bar  target  at  0.65  cycles/mrad  bar  spatial 
frequency  and  its  images  including  and  excluding  aliasing  effects  are  plotted  in  Figure 
K.4.  Contrast  enhancement  due  to  aliasing  is  noticeable. 


An  analytical  approach  to  this  question  gives  a  result  similar  to  that  obtained 
using  virtual  thermal  image  processing  model.  Equation  4  gives  the  horizontal  MRTD 
target  spectrum,  which  is  the  two-dimensional  Fourier  transform  of  the  four-bar  pattern. 


'  i  A 


\rfxj 


{-sm{nfxW)+  sin(3;zf)IT)-  sin(5/zf)IT)+ sin(7^f  W)) 


xf\W 


sm(l4sW)  (4) 


where  fx  and  fy  are  the  spatial  frequencies  along  the  horizontal  and  vertical  directions.  A 
plot  of  this  spectrum  along  the  x-axis  (fy=0)  is  given  in  Figure  K.5.  Note  the  repetition  of 
the  spectrum  at  intervals  of  the  sampling  frequency.  The  figure  shows  that  the  sample 
generated  replicas  of  the  original  spectrum  overlap  with  the  target  spectrum  (baseband) 
and  aliasing  occurs.  It  is  important  to  note  that  this  result  is  in  line  with  the  result 
obtained  from  the  virtual  thermal  image-processing  model.  Aliasing  is  an  issue  for  bar- 
targets  at  four-bar  spatial  frequencies  (1)  below  the  Nyquist  limit  (3). 
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5.  Aliasing  Effects  Above  the  Nyquist  Limit 

A  number  of  approaches  have  been  proposed  to  account  for  aliasing  in  imagery.  It 
has  been  argued  that  aliasing  can  be  adequately  represented  as  signal-dependent  additive 
noise  [Ref.  5].  Another  important  observation  made  during  the  experiments  with  the 
virtual  thermal  image-processing  model  is  that  aliasing  effects  manifest  themselves  in 
much  different  fonn  than  noise  in  imagery.  As  demonstrated  in  Figure  K.6,  effects  of 
noise  are  dominant  at  small  spatial  frequencies  and  noise  tends  to  mask  the  target  pattern. 
Aliasing  effects,  on  the  other  hand  become  more  important  at  higher  spatial  frequencies 
and  aliasing  creates  distortion  in  imagery.  The  four-bar  target  may  appear  as  a  distorted 
three-bar  pattern.  A  target  image  at  a  high  spatial  frequency  where  aliasing  is  dominant  is 
presented  in  Figure  K.7.  At  the  same  spatial  frequency,  noise  alone  does  not  distort  the 
image;  four  bars  are  resolvable  in  Figure  K.8.  These  results  suggest  that  aliasing  cannot 
be  adequately  represented  as  noise. 
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Figure  K.l.  A  Typical  MRTD  Plot  “After  [Ref.  1]”. 


Figure  K.2.  Staring  FPA  as  a  Two-Dimensional  Sampling  Lattice  “From  [Ref.  3]”. 
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Figure  K.3.  Block  Diagram  of  the  Virtual  Thermal  Image-Processing  Model. 
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Figure  K.4.  Single  Lines  Passing  through  the  Target  and  Image  Centers  (Four- 

Bar  Frequency:  0.65  cy/mrad). 


spatial  frequency  (cy/mrad) 


Figure  K.5.  Horizontal  MRTD  Bar-Pattern  along  the  x-Axis  (fy  =  0)  (Four-Bar 

Target  Frequency:  0.65  cy/mrad). 
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x  space,  mrad 


Figure  K.6.  Image  with  Noise  and  Aliasing  at  0.05  Cycles/mrad. 


x  space,  mrad 


Figure  K.7.  Image  with  Noise  and  Aliasing  at  1.10  Cycles/mrad. 
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y  space,  mrad 


x  space,  mrad 


Figure  K.8.  Image  with  only  Noise  (no  Aliasing)  at  1.10  Cycles/mrad. 
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